On the Error Estimates of a New Operator Splitting Scheme for the Navier-Stokes Equations with Coriolis Force

Abstract

An operator splitting scheme is introduced for the numerical solution of the incompressible Navier-Stokes equations with Coriolis force. Under some mild regularity assumptions on the continuous solution, error estimates and the stability analysis for the velocity and the pressure of the new operator splitting scheme are obtained. Some numerical results are presented to verify the theoretical predictions.

1. Introduction

In this paper, we consider the numerical approximation of the unsteady Navier-Stokes equations with Coriolis force: 𝜕𝐮+()𝜕𝑡𝐮⋅∇𝐮+2𝜔×𝐮−𝜈Δ𝐮+∇𝑝=𝐟,∇⋅𝐮=0, (1.1) where Ω is a bounded region in 𝑅^𝑑 (𝑑=2,3) with a sufficiently regular boundary 𝜕Ω. 𝐮 is the velocity field, 𝑝 is the pressure divided by the density (i.e., the kinematic pressure), 𝜈=1/Re is the kinematic viscosity coefficients, Re is the Reynolds number, 𝐟 is the vector of body forces, 𝜔 is the angular velocity vector, 𝐫 is the radius vector from the center of coordinates, and 2𝜔×𝐮 is the so-called Coriolis force.

For the sake of completeness, the equations should be supplemented with appropriate initial and boundary condition: 𝐮(𝐱,0)=𝑢₀(𝐱)𝐱∈Ω,𝐮(𝐱,𝑡)=𝑔(𝑡,𝐱)𝐱∈𝜕Ω,𝑡∈[0,𝑇]. (1.2) The difficulties for the numerical simulation of the incompressible flows are mainly of two kinds: nonlinearity and incompressibility. The velocity and the pressure are coupled by the incompressibility constraint, which requires that the solution spaces, to which the velocity and the pressure belong, verify the so-called inf-sup condition. To overcome these difficulties, operator splitting method, which can be viewed as the fractional step method, is introduced. Fractional step methods allow to separate the effects of the different operators appearing in the equation by splitting the time advancement into a series of substeps. In addition, the cost of simulation can be also reduced by using the fractional step method. However, these methods have a main disadvantage that splitting error is inevitable unless the operator is commute.

The origin of this category of methods is contributed to the work of Chorin and Témam [1, 2]. They developed the so-called projection method, in which the second step consists of the projection of an intermediate velocity field onto the space of solenoidal vector field. The most attractive feature of projection methods is that, at each time step, one only needs to solve a sequence of decoupled elliptic equations for the velocity and the pressure, which makes it very efficient for large-scale numerical simulations. However, several issues related to these methods still deserve further analysis, and perhaps the most salient of these are the behavior of the computed pressure near boundaries and the stability of the pressure itself. The incompatibility of the projection boundary conditions may introduce a numerical boundary layer of size 𝑂(√𝜈Δ𝑡) [3, 4], where 𝜈 is the kinematic viscosity and Δ𝑡 is the time step size. The end-of-step velocities of the projection do not converge in the space 𝐻¹₀(Ω)^𝑑, since they do not satisfy the correct boundary conditions.

These methods have been widely investigated. Guermond et al. in [5] review theoretical and numerical convergence results available for projection methods. In [6–9], the analysis on first-order accurate schemes in the time size is presented. In [10, 11], Shen derived a second-order error estimates for the projection method. Olshanskii et al. [12] proposed a projection method for the Navier-Stokes equations with Coriolis force and study the accuracy of its semidiscrete form. In [13], a new discrete projection method for the numerical solution of the Navier-Stokes equations with Coriolis force is presented, where the scheme is treated as an incomplete LU factorization of the transition operator for fully implicit time discretization. In [14], complex 3D simulations of the Stirred Tank Reactor model by a modified discrete projection method for the rotating incompressible flow are presented. Numerical experiments from [13, 14] show that including 𝜔-term in the second step enhances stability and accuracy of the scheme for the case of dominating Coriolis forces.

In this paper, we will consider the unsteady Navier-Stokes equations with Coriolis force (1.1). Using the technique developed in [7, 11] for the case of 𝜔=0, a new scheme is introduced, which is a two-step scheme and allows to enforce the original boundary conditions of the problem in all substeps of the scheme. Some error estimates of both velocity and pressure for the proposed operator splitting scheme are given, which leads to the convergence of both the intermediate and the end-of-step velocities of the method to a continuous solution in the spaces 𝐿²(Ω)^𝑑 and 𝐻¹₀(Ω)^𝑑 as in [15].

The remainder of this paper is organized as follows. In Section 2, we introduce some notations and assumptions, such as the regularity assumption for their solution. In Section 3, we present a new operator splitting scheme. In Section 4, the stability analysis is presented, then in Section 5, some error estimates for the intermediate, end-of-step velocity, and the pressure are given. Finally, in Section 6, some numerical results are presented to illustrate the theoretical results.

2. Function Setting

In order to study approximation scheme for problem (1.1), the following notations and assumptions are presented. we denote by (⋅,⋅) and ‖⋅‖ the inner product and norm in 𝐿²(Ω) or 𝐿²(Ω)^𝑑. The spaces 𝐻¹₀(Ω) and 𝐻¹₀(Ω)^𝑑 are equipped with their usual norm; that is, ‖𝐮‖²₁=_Ω|∇𝑢(𝐱)|²d𝐱. (2.1) The norm in 𝐻^𝑠(Ω) will be denoted simply by ‖⋅‖_𝑠. We will use ⟨⋅,⋅⟩ to denote the duality between 𝐻^−𝑠(Ω) and 𝐻^𝑠₀(Ω) for all 𝑠>0.

The following subspace is also introduced: 𝑉=𝐮∈𝐻¹₀(Ω)^𝑑:div𝐮=0,𝐻=𝐮∈𝐿²(Ω)^𝑑:div𝐮=0,𝐮⋅𝐧=0. (2.2) For the treatment of the convective term, the following trilinear form is given: 𝑏(𝐮,𝐯,𝐰)=_Ω(𝐮⋅∇)𝐯⋅𝐰d𝐱. (2.3) It is well known that 𝑏(⋅,⋅,⋅) is continuous in 𝐻^𝑚1(Ω)×𝐻^𝑚2+1(Ω)×𝐻^𝑚3(Ω), provided 𝑚₁+𝑚₂+𝑚₃≥𝑑/2 if 𝑚_𝑖≠𝑑/2,𝑖=1,2,3, and this form is skew-symmetric with respect to its last two arguments, that is, 𝑏(𝐮,𝐯,𝐰)=−𝑏(𝐮,𝐰,𝐯)∀𝐮∈𝐻,𝐯,𝐰∈𝐻¹₀(Ω)^𝑑. (2.4) In particular, we have 𝑏(𝐮,𝐯,𝐯)=0∀𝐮∈𝐻,𝐯∈𝐻¹₀(Ω)^𝑑, (2.5) and for 𝑑≤4, 𝑏()≤𝑐‖𝐮‖𝐮,𝐯,𝐰¹‖𝐯‖₁‖𝐰‖₁,𝑐‖𝐮‖‖𝐯‖₂‖𝐰‖₁,𝑐‖𝐮‖₁‖𝐯‖₂‖𝐰‖,𝑐‖𝐮‖‖𝐯‖₁‖𝐰‖₂,𝑐‖𝐮‖₂‖𝐯‖₁‖𝐰‖,𝑐‖𝐮‖₁‖𝐯‖₁‖𝐰‖^1/2‖𝐰‖₁^1/2,𝑐‖𝐮‖^1/2‖𝐮‖₁^1/2‖𝐯‖₁‖𝐰‖₁. (2.6) We also define the Stokes operator: 𝐴𝐮=−𝑃_𝐻Δ𝐮,∀𝐮∈𝐷(𝐴)=𝑉∩𝐻²(Ω)^𝑑, (2.7) where 𝑃_𝐻 is an orthogonal projector in the Hilbert space 𝐿²(Ω)^𝑑 onto its subspace 𝐻. The Stokes operator 𝐴 is an unbounded positive self-adjoint closed operator in 𝐻 with domain 𝐷(𝐴), and its inverse 𝐴⁻¹ is compact in 𝐻. Having the following properties: there exists constant 𝑐₁,𝑐₂>0, such that ∀𝐮∈𝐻, 𝐴⁻¹𝐮_𝑠≤𝑐₁‖𝐮‖_𝑠−2𝑐for𝑠=1,2,₂‖𝐮‖²₋₁≤𝐴⁻¹𝐮,𝐮≤𝑐²₁‖𝐮‖²₋₁. (2.8) Furthermore, from (2.8), we will use (𝐴⁻¹𝐮,𝐮)^1/2 as an equivalent norm of 𝐻⁻¹(Ω)^𝑑 for 𝐮∈𝐻.

For the purpose of this paper, we also need the following regularity assumptions:  (A1) 𝐮₀∈𝐻¹(Ω)^𝑑∩𝑉,𝐟∈𝐿^∞(0,𝑇;𝐿²(Ω)^𝑑)∩𝐿²(0,𝑇;𝐻¹(Ω)^𝑑). In the three-dimension case, we assume additionally  (A2) sup_{𝑡∈[0,𝑇]}‖𝐮(𝑡)‖₁≤𝑀, where (A2) is automatically satisfied with some appropriate constant 𝑀 when 𝑑=2.

Under the regularity assumption (A1)-(A2), one can show that [9]  (a) sup_{𝑡∈[0,𝑇]}{‖𝐮(𝑡)‖₂+‖𝐮_𝑡(𝑡)‖+‖∇𝑝(𝑡)‖}≤𝑀₁,  (b) ∫^𝑇₀‖𝐮_𝑡(𝑡)‖²₁≤𝑀₁.

In addition, if we also assume that  (A3) 𝑓_𝑡∈𝐿²(0,𝑇;𝐿²(Ω))^𝑑 holds, we have  (c) ∫^𝑇₀‖𝐮_𝑡𝑡‖²₋₁𝑑𝑡≤𝑀₁

which will be used in the sequel. Indeed, the estimates (a-b) and the estimate (c) were proved for the Navier-Stokes without Coriolis term in [9, 16], respectively. However, adding linear skew-symmetric term 𝜔×𝐮 to the momentum equation does not change arguments from [16], but leads to (a)–(c) with constant 𝑀₁ depending, in general, on Ω   [12]. Next, we cite the following lemma, which will be frequently used.

Lemma 2.1 2.1 (Discrete Gronwall Lemma). Let 𝑦^𝑛, ℎ^𝑛, 𝑔^𝑛, and 𝑓^𝑛 be nonnegative sequences satisfying 𝑦^𝑚+Δ𝑡^𝑛=𝑚_𝑛=0ℎ^𝑛≤𝐵+Δ𝑡^𝑛=𝑚_𝑛=0𝑔^𝑛𝑦^𝑛+𝑓^𝑛,withΔ𝑡^{𝑛=[𝑇/Δ𝑡]}_𝑛=0𝑔^𝑛𝑇≤𝑀,∀0≤𝑚≤.Δ𝑡 (2.9) Assume Δ𝑡𝑔^𝑛<1 and let 𝜎=max_{0≤𝑛≤[𝑇/Δ𝑡]}(1−Δ𝑡^𝑛)⁻¹, then 𝑦^𝑚+Δ𝑡^𝑛=𝑚_𝑛=0ℎ^𝑛≤exp(𝜎𝑀)𝐵+Δ𝑡^𝑛=𝑚_𝑛=0𝑓^𝑛𝑇∀𝑚≤.Δ𝑡 (2.10) Hereafter, we will use 𝑐 to denote a generic constant which depends only on Ω,𝜈,𝑇, and constants from various Sobolev inequalities. We will denote 𝑀 as a generic positive constant which may additionally depend on 𝐮₀, 𝐟, 𝜔.

3. New Operator Splitting Scheme

Equation (1.1) can be written as 𝜕𝐮𝜕𝑡=𝐴₁+𝐴₂, (3.1) such that 𝐴₁()1=−𝐮⋅∇𝐮+2𝜈Δ𝐮,𝐴₂=12𝜈Δ𝐮−∇𝑝−2𝜔×𝐮+𝐟. (3.2) So an algorithm can be formulated as follows: for 𝑡∈[𝑡_𝑛,𝑡_𝑛+1]𝜕𝐮𝜕𝑡=𝐴₁,𝐮𝑡_𝑛,𝐱=𝐮𝑡_𝑛,→𝜕𝐮,𝐱𝜕𝑡=𝐴₂,𝐮𝑡_𝑛,𝐱=𝐮𝑡_𝑛+1,𝐱, (3.3) and take 𝐮(𝑡_𝑛+1,𝐱)≈𝐮 as the approximate solution of (1.1) at time 𝑡_𝑛+1.

The scheme (3.3) has an irreducible splitting error of order 𝑂(Δ𝑡). Hence, using a higher-order time stepping scheme does not improve the overall accuracy. So a first-order accurate semidiscrete version can be obtained as follows: let 𝐮⁰=𝐮₀, we solve successively 𝐮^𝑛+1 and {𝐮^𝑛+1,𝑝^𝑛+1} by 𝐮^𝑛+1−𝐮^𝑛+𝐮Δ𝑡^𝑛⋅∇𝐮^𝑛+1−12𝜈Δ𝐮^𝑛+1𝐮=0,_𝑛+1|_𝜕Ω=𝐠𝑡_𝑛+1,𝐮,𝐱_𝑛+1−𝐮^𝑛+1−1Δ𝑡2𝜈Δ𝐮^𝑛+1+∇𝑝^𝑛+1+2𝜔×𝐮^𝑛+1=𝐟𝑡_𝑛+1,div𝐮^𝑛+1𝐮=0,_𝑛+1=𝐠𝑡_𝑛+1,𝐱. (3.4)(3.5) Note that we have omitted the dependency to 𝐱 of the function 𝐟 to simplify our notations; we will do so for {𝐮,𝑝}.

As can be seen in (3.5), the main difference between this method and the standard projection method is the introduction of a viscous term in the incompressible step, which allows the imposition of the original boundary condition (2.6) on the end-of-step velocity 𝐮^𝑛+1. Similar ideas can be found in the 𝜃-scheme [17] in which viscosity and incompressibility are also coupled, and some other methods such as [18–20], all of which involve an incompressible step with part of the viscous term. It leads to convergence of both the intermediate and the end-of-step velocities of the method to a continuous solution in space 𝐿²(Ω) and 𝐻¹(Ω). In comparison with the 𝜃 scheme, our scheme is two steps instead of three steps in 𝜃 scheme. Moreover, the fact that 𝐮^𝑛+1 satisfies the correct boundary conditions will allow us to obtain improved error estimates comparative with the standard projection method.

Denoting the corresponding right-hand side by 𝐟, at each time step, we have to solve the following two subproblems: ()𝐮𝛼𝐮+𝐰⋅∇𝐮+𝛽Δ𝐮=𝐟,_𝑛+1|_𝜕Ω=𝐠𝑡_𝑛+1,𝐮,𝐱𝛾𝐮+𝛽Δ𝐮+∇𝑝=𝐟,div^𝑛+1𝐮=0,_𝑛+1=𝐠𝑡_𝑛+1,𝐱, (3.6)(3.7) with 𝛼=1/Δ𝑡,𝛽=−(1/2)𝜈,𝛾=1/Δ𝑡+2×𝜔.

The first step of the method is a linearized elliptic problem, which can be seen as a linear Burger’s problem. The second step is a generalized Stokes problem. To solve problem (3.6), the fixed point iterative technique is used as in [21, 22], which is cheaper than the conjugate gradient method used by the least-square technique in the corresponding advective subproblems appearing in Glowinski’s 𝜃-scheme. To solve problem (3.7), the efficient technique of the functional equation satisfied by the pressure is used, which is used in the Glowinski’s 𝜃-scheme for the corresponding Stokes problem; that is, conjugate gradient method is applied on the variational formulation of such an equation. One defect of this method is that the discrete inf-sup compatibility condition should be satisfied.

4. Stability Analysis

For the sake of simplicity, we will only consider the homogeneous boundary condition 𝐮(𝐱,𝑡)|_𝜕Ω=0, that is, 𝐠(𝐭,𝐱)=0 for the scheme (3.4)-(3.5).

Theorem 4.1. Under the assumptions (A1)-(A2), there exists a constant 𝑐, such that 𝐮^𝑛∈𝐿²0,𝑇,𝐻¹(Ω)^𝑑∩𝐿^∞0,𝑇,𝐿²(Ω)^𝑑,𝐮^𝑛∈𝐿²0,𝑇,𝐻¹(Ω)^𝑑∩𝐿^∞0,𝑇,𝐿²(Ω)^𝑑. (4.1)

Proof. We take the inner product of (3.4) with 2Δ𝑡𝐮^𝑛+1 to get 𝐮^𝑛+1−𝐮^𝑛2+𝐮^𝑛+12−𝐮^𝑛2+𝜈Δ𝑡𝐮^𝑛+12₁=0. (4.2) Next, taking the inner product of (3.5) with 2Δ𝑡𝐮^𝑛+1, and using the condition div𝐮^𝑛+1=0, we obtain 𝐮^𝑛+1−𝐮^𝑛+12+𝐮^𝑛+12−𝐮^𝑛+12+𝜈Δ𝑡𝐮^𝑛+12₁=2Δ𝑡𝐟𝑡_𝑛+1,𝐮^𝑛+1. (4.3) Combing (4.2) with (4.3), and using the Young’s inequality, we have 𝐮^𝑛+12+𝐮^𝑛+1−𝐮^𝑛+12+𝐮^𝑛+1−𝐮^𝑛2+𝜈Δ𝑡𝐮^𝑛+12₁+𝐮^𝑛+12₁≤𝑐Δ𝑡(𝐟𝑡_𝑛+1²+Δ𝑡𝐮^𝑛+12+𝐮^𝑛2. (4.4) Summing up the inequality (4.4) for 𝑛=0,…,𝑟≤𝑁, we obtain 𝐮^𝑟+12+^𝑛=𝑟_𝑛=0𝐮^𝑛+1−𝐮^𝑛+12+^𝑛=𝑟_𝑛=0𝐮^𝑛+1−𝐮^𝑛2+Δ𝑡𝜈^𝑛=𝑟_𝑛=0𝐮^𝑛+12₁+𝜈Δ𝑡^𝑛=𝑟_𝑛=0𝐮^𝑛+12₁≤𝑐Δ𝑡^𝑛=𝑟_𝑛=0𝐟(𝑡_𝑛+1)²₀+𝐮₀²+Δ𝑡^𝑛=𝑟_𝑛=0𝐮^𝑛+12. (4.5) Applying the discrete Gronwall lemma to the above inequality, we obtain 𝐮^𝑟+12+^𝑛=𝑟_𝑛=0𝐮^𝑛+1−𝐮^𝑛+12+^𝑛=𝑟_𝑛=0𝐮^𝑛+1−𝐮^𝑛2+Δ𝑡𝜈^𝑛=𝑟_𝑛=0𝐮^𝑛+12₁+𝜈Δ𝑡^𝑛=𝑟_𝑛=0𝐮^𝑛+12₁≤𝑐Δ𝑡^𝑛=𝑟_𝑛=0𝐟𝑡_𝑛+1²₀+𝐮₀². (4.6) Thus, using the regularity properties of the continuous solution 𝐮, for arbitrary 𝑛, we have 𝐮^𝑛+12+^𝑛=𝑁_𝑛=0𝐮^𝑛+1−𝐮^𝑛+12+𝐮^𝑛+1−𝐮^𝑛2+𝜈Δ𝑡^𝑛=𝑁_𝑛=0𝐮^𝑛+12₁+𝐮^𝑛+12₁≤𝑐, (4.7) which means that 𝐮^𝑛∈𝐿²0,𝑇,𝐻¹(Ω)^𝑑,𝐮^𝑛∈𝐿²0,𝑇,𝐻¹(Ω)^𝑑∩𝐿^∞0,𝑇,𝐿²(Ω)^𝑑. (4.8) From (4.2), it yields𝐮^𝑛+12≤𝐮^𝑛2. (4.9) So we have 𝐮^𝑛∈𝐿^∞(0,𝑇,𝐿²(Ω)^𝑑). The proof is completed.

Remark 4.2. The formula (4.7) can be viewed as the discrete version of the classical energy estimate for the Navier-Stokes equations [23]. From (4.6), we have 𝐮^𝑁+12≤𝑐Δ𝑡^𝑛=𝑁_𝑛=0𝑓𝑡_𝑛+1²₀+𝐮₀². (4.10) This estimate provides a meaningful bound for ‖𝐮^𝑁+1‖² for the first few time steps, that is, for small 𝑇.

5. Error Analysis

In this section, we present an error analysis of the operator splitting scheme introduced in the previous section. Firstly, we define the semidiscrete velocity error as 𝐞̂^𝑛+1=𝐮𝑡_𝑛+1−𝐮^𝑛+1,𝐞̃^𝑛+1=𝐮𝑡_𝑛+1−𝐮^𝑛+1. (5.1) We give a first estimate for 𝐞̂^𝑛+1 and 𝐞̃^𝑛+1 which shows that 𝐮^𝑛+1 and 𝐮^𝑛+1 are both order 1/2 approximations to 𝐮 in 𝑙^∞(𝐿²(Ω)^𝑑) and 𝑙²(𝐻¹(Ω)^𝑑).

Theorem 5.1. Under the regularity assumptions (A1)–(A3), there exists a constant 𝑀, such that 𝐞̂^𝑁+12+𝐞̃^𝑁+12+^𝑛=𝑁_𝑛=0𝐞̃^𝑛+1−𝐞̂^𝑛2+𝐞̂^𝑛+1−𝐞̃^𝑛+12+𝜈Δ𝑡𝐞̂^𝑛+12₁+𝐞̃^𝑛+12₁≤𝑀Δ𝑡. (5.2)

Proof. Let 𝑅^𝑛 be the truncation error defined by 1𝐮𝑡Δ𝑡^𝑛+1−𝐮𝑡_𝑛−𝜈Δ𝐮𝑡_𝑛+1+2𝜔×𝐮𝑡_𝑛+1+𝐮𝑡_𝑛+1⋅∇𝐮𝑡_𝑛+1+∇𝑝𝑡_𝑛+1=𝐟𝑡_𝑛+1+𝐑^𝑛, (5.3) where 𝐑^𝑛 is the integral residual of the Taylor series, that is, 𝐑^𝑛=1Δ𝑡^𝑡𝑛+1_𝑡𝑛𝑡−𝑡_𝑛𝐮_𝑡𝑡𝑑𝑡. (5.4) By subtracting (3.4) from (5.3), we obtain 𝐞̃^𝑛+1−𝐞̂^𝑛−1Δ𝑡2𝜈Δ𝐞̃^𝑛+1−12𝜈Δ𝐮𝑡_𝑛+1+2𝜔×𝐮𝑡_𝑛+1=𝐮^𝑛⋅∇𝐮^𝑛+1−𝐮𝑡_𝑛+1⋅∇𝐮𝑡_𝑛+1+𝐑^𝑛−∇𝑝𝑡_𝑛+1+𝐟𝑡_𝑛+1. (5.5) The nonlinear terms on the right-side can be split into three terms: 𝐮^𝑛⋅∇𝐮^𝑛+1−𝐮𝑡_𝑛+1⋅∇𝐮𝑡_𝑛+1=−𝐞̂^𝑛⋅∇𝐮^𝑛+1+𝐮𝑡_𝑛−𝐮𝑡_𝑛+1⋅∇𝐮^𝑛+1−𝐮𝑡_𝑛+1⋅∇𝐞̃^𝑛+1. (5.6) Taking the inner product of (5.5) with 2Δ𝑡𝐞̃^𝑛+1, using the identity (𝑎−𝑏,2𝑎) =|𝑎|²+|𝑎+𝑏|²−|𝑏|², we obtain 𝐞̃^𝑛+12+𝐞̃^𝑛+1−𝐞̂^𝑛2−𝐞̂^𝑛2+𝜈Δ𝑡𝐞̃^𝑛+12₁−𝜈Δ𝑡Δ𝐮𝑡_𝑛+1,𝐞̃^𝑛+1+4Δ𝑡𝜔×𝐮𝑡_𝑛+1,𝐞̃^𝑛+1=2Δ𝑡<𝐑^𝑛,𝐞̃^𝑛+1>−2Δ𝑡∇𝑝𝑡_𝑛+1,𝐞̃^𝑛+1−2Δ𝑡𝑏𝐞̂^𝑛,𝐮^𝑛+1,𝑒̃^𝑛+1+2Δ𝑡𝑏𝐮𝑡_𝑛−𝐮𝑡_𝑛+1,𝐮^𝑛+1,𝐞̃^𝑛+1−2Δ𝑡𝑏𝐮𝑡_𝑛+1,𝐞̃^𝑛+1,𝐞̃^𝑛+1+2Δ𝑡𝐟𝑡_𝑛+1,𝐞̃^𝑛+1. (5.7) On the other hand, we derive from (3.5) that 𝐞̂^𝑛+1−𝐞̃^𝑛+1+1Δ𝑡2𝜈Δ𝐮^𝑛+1−∇𝑝^𝑛+1−2𝜔×𝐮^𝑛+1+𝐟𝑡_𝑛+1=0. (5.8) Taking the inner product of the last equality with 2Δ𝑡𝐞̂^𝑛+1, and using div𝐞̂^𝑛+1=0, we obtain 𝐞̂^𝑛+12+𝐞̂^𝑛+1−𝐞̃^𝑛+12−𝐞̃^𝑛+12+Δ𝑡𝜈Δ𝐮^𝑛+1,𝐞̂^𝑛+1−4Δ𝑡𝜔×𝐮^𝑛+1,𝐞̂^𝑛+1+2Δ𝑡𝐟𝑡_𝑛+1,𝐞̂^𝑛+1=0. (5.9) Combing (5.7) with (5.9), we obtain 𝐞̂^𝑛+12+𝐞̃^𝑛+1−𝐞̂^𝑛2+𝐞̂^𝑛+1−𝐞̃^𝑛+12+Δ𝑡𝜈𝐞̂^𝑛+12₁+𝜈Δ𝑡𝐞̃^𝑛+12₁=2Δ𝑡<𝐑^𝑛,𝐞̃^𝑛+1>−2Δ𝑡∇𝑝𝑡_𝑛+1,𝐞̃^𝑛+1−2Δ𝑡𝑏𝐞̂^𝑛,𝐮^𝑛+1,𝐞̃^𝑛+1+2Δ𝑡𝑏𝐮𝑡_𝑛−𝐮𝑡_𝑛+1,𝐮^𝑛+1,𝐞̃^𝑛+1−2Δ𝑡𝑏𝐮𝑡_𝑛+1,𝐞̃^𝑛+1,𝐞̃^𝑛+1+2Δ𝑡𝐟𝑡_𝑛+1,𝐞̃^𝑛+1−𝐞̂^𝑛+1+Δ𝑡𝜈Δ𝐮𝑡_𝑛+1,𝐞̃^𝑛+1−𝐞̂^𝑛+1+𝐞̂^𝑛2−4Δ𝑡𝜔×𝐮𝑡_𝑛+1,𝐞̃^𝑛+1−𝐞̂^𝑛+1. (5.10) We bound each term in the right-hand side of (5.10) independently.
Taylor residual term: 2Δ𝑡𝐑^𝑛,𝐞̃^𝑛+1≤2Δ𝑡𝐑^𝑛₋₁𝐞̃^𝑛+1₁≤Δ𝑡𝜈6𝐞̃^𝑛+12₁+𝑐Δ𝑡^2𝑡𝑛+1_𝑡𝑛𝐮_𝑡𝑡²₋₁𝑑𝑡. (5.11) Pressure gradient term: −2Δ𝑡∇𝑝𝑡_𝑛+1,𝐞̃^𝑛+1=−2Δ𝑡∇𝑝𝑡_𝑛+1,𝐞̃^𝑛+1−𝐞̂^𝑛≤12𝐞̃^𝑛+1−𝐞̂^𝑛2+𝑐Δ𝑡²∇𝑝(𝑡_𝑛+1)². (5.12) Nonlinear term: −2Δ𝑡𝑏𝐞̂^𝑛,𝐮^𝑛+1,𝐞̃^𝑛+1=−2Δ𝑡𝑏𝐞̂^𝑛,𝐮𝑡_𝑛+1,𝐞̃^𝑛+1≤𝑐Δ𝑡𝐞̂^𝑛𝐞̃^𝑛+1₁𝐮(𝑡_𝑛+1)₂≤Δ𝑡𝜈6𝐞̃^𝑛+12₁+𝑐Δ𝑡𝐞̂^𝑛2,2Δ𝑡𝑏𝐮𝑡_𝑛−𝐮𝑡_𝑛+1,𝐮^𝑛+1,𝐞̃^𝑛+1=−2Δ𝑡𝑏𝐮𝑡_𝑛−𝐮𝑡_𝑛+1,𝐞̃^𝑛+1,𝐮^𝑛+1=−2Δ𝑡𝑏𝐮𝑡_𝑛−𝐮𝑡_𝑛+1,𝐞̃^𝑛+1,𝐮𝑡_𝑛+1≤𝑐Δ𝑡𝐮(𝑡_𝑛)−𝐮(𝑡_𝑛+1)𝐞̃^𝑛+1₁𝐮(𝑡_𝑛+1)₂≤Δ𝑡𝜈6𝐞̃^𝑛+12₁+𝑐Δ𝑡^2𝑡𝑛+1_𝑡𝑛𝐮_𝑡²𝐮𝑡𝑑𝑡,−2Δ𝑡𝑏_𝑛+1,𝐞̃^𝑛+1,𝐞̃^𝑛+1=0. (5.13) The external term: 2Δ𝑡𝐟𝑡_𝑛+1,𝐞̃^𝑛+1−𝐞̂^𝑛+1≤16𝐞̃^𝑛+1−𝐞̂^𝑛+12+𝑐Δ𝑡²𝐟𝑡_𝑛+1². (5.14) The viscous term: 2Δ𝑡𝜈Δ𝐮𝑡_𝑛+1,𝐞̃^𝑛+1−𝐞̂^𝑛+1≤16𝐞̃^𝑛+1−𝐞̂^𝑛+12+𝑐Δ𝑡²Δ𝐮𝑡_𝑛+1²≤16𝐞̃^𝑛+1−𝐞̂^𝑛+12+𝑐Δ𝑡²𝐮𝑡_𝑛+1²₂. (5.15) The rotating term: −4Δ𝑡𝜔×𝐮𝑡_𝑛+1,𝐞̃^𝑛+1−𝐞̂^𝑛+1≤𝑐Δ𝑡²𝜔×𝐮𝑡_𝑛+1²+16𝐞̃^𝑛+1−𝐞̂^𝑛+12≤𝑐Δ𝑡²‖𝜔‖²𝐮𝑡_𝑛+1²+16𝐞̃^𝑛+1−𝐞̂^𝑛+12. (5.16) Inserting the above estimates into (5.10), we obtain 𝐞̂^𝑛+12+12𝐞̃^𝑛+1−𝐞̂^𝑛2+12𝐞̂^𝑛+1−𝐞̃^𝑛+12+𝜈Δ𝑡𝐞̂^𝑛+12₁+12𝜈Δ𝑡𝐞̃^𝑛+12₁≤𝑐Δ𝑡^2𝑡𝑛+1_𝑡𝑛𝐮_𝑡𝑡²₋₁𝑑𝑡+𝑐Δ𝑡²∇𝑝𝑡_𝑛+1²+𝑐Δ𝑡𝐞̂^𝑛2+𝑐Δ𝑡^2𝑡𝑛+1_𝑡𝑛𝐮_𝑡²𝑑𝑡+𝑐Δ𝑡²𝐟𝑡_𝑛+1²+𝑐Δ𝑡²𝐮𝑡_𝑛+1²₂+𝑐Δ𝑡²‖𝜔‖²𝐮𝑡_𝑛+1². (5.17) Summing up the inequality (5.17) for 𝑛=0,…,𝑁, we get 𝐞̂^𝑁+12+^𝑛=𝑁_𝑛=012𝐞̃^𝑛+1−𝐞̂^𝑛2+12𝐞̂^𝑛+1−𝐞̃^𝑛+12+𝜈Δ𝑡𝐞̂^𝑛+12₁+12𝜈Δ𝑡𝐞̃^𝑛+12₁≤𝑐Δ𝑡Δ𝑡^𝑡𝑇₀𝐮_𝑡𝑡²₋₁𝑑𝑡+sup_{𝑡∈[0,𝑇]}∇𝑝(𝑡_𝑛+1)²+Δ𝑡^𝑇₀𝐮_𝑡²+𝑑𝑡_{𝑛=𝑁𝑛=0}Δ𝑡𝐟𝑡_𝑛+1²+^𝑛=𝑁_𝑛=0Δ𝑡𝐮𝑡_𝑛+12+^𝑛=𝑁_𝑛=0Δ𝑡‖𝜔‖²𝐮𝑡_𝑛+1². (5.18) By applying the discrete Gronwall lemma to the above inequality, we derive 𝐞̂^𝑁+12+^𝑛=𝑁_𝑛=012𝐞̃^𝑛+1−𝐞̂^𝑛2+12𝐞̂^𝑛+1−𝐞̃^𝑛+12+𝜈Δ𝑡𝐞̂^𝑛+12₁+12𝜈Δ𝑡𝐞̃^𝑛+12₁≤𝑐Δ𝑡Δ𝑡^𝑇₀𝐮_𝑡𝑡²₋₁𝑑𝑡+sup_{[]𝑡∈0,𝑇}∇𝑝𝑡_𝑛+1²+Δ𝑡^𝑇₀𝐮_𝑡²+𝑑𝑡_{𝑛=𝑁𝑛=0}Δ𝑡𝐟𝑡_𝑛+1²+^𝑛=𝑁_𝑛=0Δ𝑡𝐮𝑡_𝑛+1²₂+^𝑛=𝑁_𝑛=0Δ𝑡‖𝜔‖²𝐮𝑡_𝑛+1². (5.19) Using the regularity properties of 𝐮, we obtain 𝐞̂^𝑁+12+^𝑛=𝑁_𝑛=0𝐞̃^𝑛+1−𝐞̂^𝑛2+𝐞̂^𝑛+1−𝐞̃^𝑛+12+𝜈Δ𝑡𝐞̂^𝑛+12₁+𝜈Δ𝑡𝐞̃^𝑛+12₁≤𝑀Δ𝑡. (5.20) Finally, the bounds for 𝐮^𝑛+1 follow from (5.20) and the triangle inequality. Theorem 5.1 is proved.

Remark 5.2. Theorem 5.1 shows that 𝐮^𝑛+1₁≤𝑀,𝐮^𝑛+1₁≤𝑀, (5.21) since ‖𝐞̂^𝑛+1‖₁≤𝑀,‖𝐞̃^𝑛+1‖₁≤𝑀. Moreover, we also have 𝐮^𝑛+1₀≤𝑀Δ𝑡^1/2,𝐮^𝑛+1₀≤𝑀Δ𝑡^1/2. (5.22)
Next, we will use the previous result to improve the error estimates for the velocity and give an error estimate for the pressure as well. The result shows that both 𝐮^𝑛+1 and 𝐮^𝑛+1 are weakly first-order approximations to 𝐮(𝑡_𝑛+1) in 𝑙²(𝐿²(Ω)^𝑑).

Theorem 5.3. Under the regularity assumptions (A1)–(A3), there exists a constant 𝑀, such that 𝐞̂^𝑁+12₋₁+^𝑁_𝑁=0𝐞̂^𝑛+1−𝐞̂^𝑛2₋₁+𝜈Δ𝑡𝑒̃^𝑛+12+𝜈Δ𝑡𝑒̂^𝑛+12≤𝑀Δ𝑡²,Δ𝑡^𝑛=𝑁_𝑛=0𝑝(𝑡_𝑛+1)−𝑝^𝑛+1_𝐿20Ω≤𝑀Δ𝑡. (5.23)(5.24)

Proof. Taking the sum of (3.4) and (3.5), we obtain 𝐮^𝑛+1−𝐮^𝑛+𝐮Δ𝑡^𝑛⋅∇𝐮^𝑛+1−12𝜈Δ𝐮^𝑛+1−12𝜈Δ𝐮^𝑛+1+∇𝑝^𝑛+1+2𝜔×𝐮^𝑛+1=𝐟𝑡_𝑛+1. (5.25) Let us denote 𝑞^𝑛+1=𝑝𝑡_𝑛+1−𝑝^𝑛+1. (5.26) Subtracting (5.25) from (5.3), we obtain 𝐞̂^𝑛+1−𝐞̂^𝑛−1Δ𝑡2𝜈Δ𝐞̃^𝑛+1−12𝜈Δ𝐞̂^𝑛+1+∇𝑞^𝑛+1+2𝜔×𝐞̂^𝑛+1=𝐮^𝑛⋅∇𝐮^𝑛+1−𝐮𝑡_𝑛+1⋅∇𝐮𝑡_𝑛+1+𝐑^𝑛. (5.27)Taking the inner product of the expression (5.27) with 2Δ𝑡𝐴⁻¹𝐞̂^𝑛+1, we obtain 𝐞̂^𝑛+1,𝐴⁻¹𝐞̂^𝑛+1−𝐞̂^𝑛,𝐴⁻¹𝐞̂^𝑛+𝐞̂^𝑛+1−𝐞̂^𝑛,𝐴⁻¹𝐞̂^𝑛+1−𝐞̂^𝑛−𝜈Δ𝑡Δ𝐞̃^𝑛+1,𝐴⁻¹𝐞̂^𝑛+1−𝜈Δ𝑡Δ𝐞̂^𝑛+1,𝐴⁻¹𝐞̂^𝑛+1+4Δ𝑡𝜔×𝑒̂^𝑛+1,𝐴⁻¹𝐞̂^𝑛+1=2Δ𝑡𝑏𝐮^𝑛,𝐮^𝑛+1,𝐴⁻¹𝐞̂^𝑛+1−2Δ𝑡𝑏𝐮𝑡_𝑛+1,𝐮𝑡_𝑛+1,𝐴⁻¹𝐞̂^𝑛+1+2Δ𝑡𝐑^𝑛,𝐴⁻¹𝐞̂^𝑛+1. (5.28) Since −𝜈Δ𝑡Δ𝐞̃^𝑛+1,𝐴⁻¹𝐞̂^𝑛+11≥𝜈Δ𝑡4𝐞̂^𝑛+12−𝑐𝐞̃^𝑛+1−𝐞̂^𝑛+12−𝜈Δ𝑡Δ𝐞̂^𝑛+1,𝐴⁻¹𝐞̂^𝑛+1=𝜈Δ𝑡𝐞̂^𝑛+12, (5.29) together with (5.28), it yields 𝐞̂^𝑛+12₋₁−𝐞̂^𝑛2₋₁+𝐞̂^𝑛+1−𝐞̂^𝑛2₋₁+54𝜈Δ𝑡𝐞̂^𝑛+12≤2Δ𝑡𝑏𝐮^𝑛,𝐮^𝑛+1,𝐴⁻¹𝐞̂^𝑛+1−2Δ𝑡𝑏𝐮𝑡_𝑛+1,𝐮𝑡_𝑛+1,𝐴⁻¹𝐞̂^𝑛+1+2Δ𝑡𝐑^𝑛,𝐴⁻¹𝐞̂^𝑛+1+𝑐Δ𝑡𝐞̃^𝑛+1−𝐞̂^𝑛+12−4Δ𝑡𝜔×𝐞̂^𝑛+1,𝐴⁻¹𝐞̂^𝑛+1. (5.30) Similar to (5.6) for the nonlinear term, together with (5.30), it yields 𝐞̂^𝑛+12₋₁−𝐞̂^𝑛2₋₁+𝐞̂^𝑛+1−𝐞̂^𝑛2₋₁+54𝜈Δ𝑡𝐞̂^𝑛+12=−2Δ𝑡𝑏𝐞̂^𝑛,𝐮^𝑛+1,𝐴⁻¹𝐞̂^𝑛+1+2Δ𝑡𝑏𝐮𝑡_𝑛−𝐮𝑡_𝑛+1,𝐮^𝑛+1,𝐴⁻¹𝐞̂^𝑛+1−2Δ𝑡𝑏𝐮𝑡_𝑛+1,𝐞̃^𝑛+1,𝐴⁻¹𝐞̂^𝑛+1+2Δ𝑡𝐑^𝑛,𝐴⁻¹𝐞̂^𝑛+1+𝑐Δ𝑡𝐞̃^𝑛+1−𝐞̂^𝑛+12−4Δ𝑡𝜔×𝐞̂^𝑛+1,𝐴⁻¹𝐞̂^𝑛+1. (5.31) We will focus on the right-hand side as follows.
The Coriolis term is estimated as follows: 4Δ𝑡𝜔×𝐞̂^𝑛+1,𝐴⁻¹𝐞̂^𝑛+1≤𝑐Δ𝑡𝐴⁻¹𝐞̂^𝑛+1₁𝐞̂^𝑛+1≤𝑐Δ𝑡𝐞̂^𝑛+1₋₁𝐞̂^𝑛+1≤𝜈Δ𝑡𝐞̂20^𝑛+12+𝑐Δ𝑡𝐞̂^𝑛+12₋₁. (5.32) For the Taylor residual term, we have 2Δ𝑡𝐑^𝑛,𝐴⁻¹𝐞̂^𝑛+1≤𝑐Δ𝑡𝐑^𝑛₋₁𝐴⁻¹𝐞̂^𝑛+1₁≤Δ𝑡𝐞̂^𝑛+12₋₁+𝑐Δ𝑡^2𝑡𝑛+1_𝑡𝑛𝐮_𝑡𝑡²₋₁𝑑𝑡. (5.33) For the nonlinear term, it yields −2Δ𝑡𝑏𝐞̂^𝑛,𝐮^𝑛+1,𝐴⁻¹𝐞̂^𝑛+1=2Δ𝑡𝑏𝐞̂^𝑛,𝐴⁻¹𝐞̂^𝑛+1,𝐮^𝑛+1=2Δ𝑡𝑏𝐞̂^𝑛,𝐴⁻¹𝐞̂^𝑛+1,𝐮𝑡_𝑛+1−2Δ𝑡𝑏𝐞̂^𝑛,𝐴⁻¹𝐞̂^𝑛+1,𝐞̃^𝑛+1=𝑇₁+𝑇₂. (5.34) Using (2.6)-(2.8), we derive 𝑇₁≤𝑐Δ𝑡𝐞̂^𝑛𝐴⁻¹𝐞̂^𝑛+1₁𝐮𝑡_𝑛+12≤𝑐Δ𝑡𝐞̂^𝑛𝐞̂^𝑛+1₋₁≤𝑐Δ𝑡𝐞̂^𝑛+1+𝐞̂^𝑛+1−𝐞̃^𝑛+1+𝐞̃^𝑛+1−𝐞̂^𝑛𝐞̂^𝑛+1₋₁≤𝜈Δ𝑡𝐞̂20^𝑛+12+𝑐Δ𝑡𝐞̂^𝑛+1−𝐞̃^𝑛+12+𝐞̃^𝑛+1−𝐞̂^𝑛2+𝑐Δ𝑡𝐞̂^𝑛+12₋₁, (5.35) and together with Theorem 5.1, 𝑇₂≤𝑐Δ𝑡𝐞̂^𝑛𝐴⁻¹𝐞̂^𝑛+1₂𝐞̃^𝑛+1₁≤𝑐Δ𝑡𝐞̂^𝑛𝐞̂^𝑛+1𝐞̃^𝑛+1₁≤𝑐Δ𝑡^3/2𝐞̂^𝑛+1𝐞̃^𝑛+1₁≤𝜈Δ𝑡𝐞̂20^𝑛+12+𝑐Δ𝑡²𝐞̃^𝑛+1₁,−2Δ𝑡𝑏𝐮𝑡_𝑛+1,𝑒̃^𝑛+1,𝐴⁻¹𝐞̂^𝑛+1=2Δ𝑡𝑏𝑢𝑡_𝑛+1,𝐴⁻¹𝐞̂^𝑛+1,𝑒̃^𝑛+1≤𝑐Δ𝑡𝐮𝑡_𝑛+12𝐴⁻¹𝐞̂^𝑛+1₁𝑒̃^𝑛+1≤𝑐Δ𝑡𝐞̂^𝑛+1₋₁𝑒̃^𝑛+1≤𝑐Δ𝑡𝐞̂^𝑛+1₋₁𝑒̂^𝑛+1−𝑒̃^𝑛+1+𝑒̂^𝑛+1≤𝜈Δ𝑡𝐞̂20^𝑛+12+𝑐Δ𝑡𝐞̂^𝑛+1−𝐞̃^𝑛+12+𝑐Δ𝑡𝐞̂^𝑛+12₋₁. (5.36) Similarly, 2Δ𝑡𝑏𝐮𝑡_𝑛−𝑢𝑡_𝑛+1,𝐮^𝑛+1,𝐴⁻¹𝐞̂^𝑛+1≤𝐮𝑡_𝑛−𝑢𝑡_𝑛+1𝐮^𝑛+1₁𝐴⁻¹𝐞̂^𝑛+1₂≤𝜈Δ𝑡𝐞̂16^𝑛+12+𝑐Δ𝑡^2𝑡𝑛+1_𝑡𝑛𝐮_𝑡²𝑑𝑡. (5.37) Inserting the above inequality into (5.31), we obtain 𝐞̂^𝑛+12₋₁−𝐞̂^𝑛2₋₁+𝐞̂^𝑛+1−𝐞̂^𝑛2₋₁+𝜈Δ𝑡𝑒̂^𝑛+12=𝑐Δ𝑡^2𝑡𝑛+1_𝑡𝑛𝐮_𝑡𝑡²₋₁𝑑𝑡+^𝑡𝑛+1_𝑡𝑛𝐮_𝑡²𝑑𝑡+𝑐Δ𝑡𝐞̂^𝑛+1−𝐞̃^𝑛+12+𝐞̃^𝑛+1−𝐞̂^𝑛2+𝑐Δ𝑡²𝐞̃^𝑛+1₁+𝑐Δ𝑡𝐞̂^𝑛+12₋₁. (5.38) Taking the sum of the above inequality for 𝑛=0,…,𝑁, using the regularity assumption of the solution 𝐮 and Theorem 5.1, yields 𝐞̂^𝑁+12₋₁+^𝑁_𝑁=0𝐞̂^𝑛+1−𝐞̂^𝑛2₋₁+𝜈Δ𝑡𝑒̂^𝑛+12≤𝑚Δ𝑡²+𝑐Δ𝑡^𝑁_𝑁=0𝐞̂^𝑛+12₋₁. (5.39) By applying the discrete Gronwall lemma to the last inequality, we obtain 𝐞̂^𝑁+12₋₁+^𝑁_𝑁=0𝐞̂^𝑛+1−𝐞̂^𝑛2₋₁+𝜈Δ𝑡𝑒̂^𝑛+12≤𝑀Δ𝑡². (5.40) For 𝐮^𝑛+1, we have 𝜈Δ𝑡^𝑁_𝑁=0𝐞̃^𝑛+12≤𝜈Δ𝑡^𝑁_𝑁=0𝐞̂^𝑛+1+𝐞̃^𝑛+1−𝐞̂^𝑛+1≤𝑀Δ𝑡², (5.41) and together with (5.40), we derive (5.23).
Next, we derive the estimate for the pressure; we recast (5.27) as ∇𝑞^𝑛+1=12𝜈Δ𝐞̃^𝑛+1+12𝜈Δ𝐞̂^𝑛+1−𝐞̂^𝑛+1−𝐞̂^𝑛+𝐮Δ𝑡_𝑛⋅∇𝐮^𝑛+1−𝐮𝑡_𝑛+1⋅∇𝐮𝑡_𝑛+1+𝐑^𝑛−2𝜔×𝐞̂^𝑛+1, (5.42) firstly, by using (2.6) and Theorem 5.3, for all 𝐯∈𝐻¹₀(Ω)^𝑑, 𝐮^𝑛⋅∇𝐮^𝑛+1−𝐮𝑡_𝑛+1⋅∇𝐮𝑡_𝑛+1≤𝑐𝐮𝑡_𝑛+1−𝐮𝑡_𝑛𝐮𝑡_𝑛+12‖𝐯‖₁+𝑐𝐞̂^𝑛₁𝐮𝑡_𝑛+11‖𝐯‖₁+𝑐𝐮^𝑛₁𝐞̃^𝑛+1₁‖𝐯‖₁≤𝐞̃^𝑛+1₁+𝐞̂^𝑛₁+𝐮𝑡_𝑛+1−𝐮𝑡_𝑛‖𝐯‖₁. (5.43) Using the Schwarz inequality, we have also, for all 𝐯∈𝐻¹₀(Ω)^𝑑, 12𝜈Δ𝐞̃^𝑛+1+12𝜈Δ𝐞̂^𝑛+1−𝐞̂^𝑛+1−𝐞̂^𝑛𝐞̂Δ𝑡−2𝜔×^𝑛+1+𝐑^𝑛≤1,𝐯𝐞̂Δ𝑡^𝑛+1−𝐞̂^𝑛₋₁+𝐑^𝑛₋₁+12𝜈𝐞̂^𝑛+1₁+12𝜈𝐞̃^𝑛+1₁+𝑐𝐞̂^𝑛+1₁‖𝐯‖₁. (5.44) Finally, we derive 𝑛𝑞^𝑛+1_𝐿20Ω≤𝑐sup_{𝑣∈𝐻10(Ω)𝑑}∇𝑞^𝑛+1,𝐯‖𝐯‖₁≤𝑐𝐞̂Δ𝑡^𝑛+1−𝐞̂^𝑛₋₁+𝑐𝐑^𝑛₋₁+𝐞̃^𝑛+1₁+𝐞̂^𝑛+1₁+𝐞̂^𝑛₁+𝐮𝑡_𝑛+1−𝐮𝑡_𝑛. (5.45) Therefore, by using Theorem 5.1 and (5.40), we obtain (5.24).

The error estimate of Theorem 5.3 can be improved to first order on the norms of 𝑙^∞(𝐿²(Ω)^𝑑) and 𝑙²(𝐻¹₀(Ω)^𝑑) for the end-of-step velocities 𝐮^𝑛+1.

Theorem 5.4. Under the regularity assumptions (A1)–(A3), for small enough Δ𝑡, there exists a constant 𝑀, such that 𝐞̂^𝑁+12+^𝑛=𝑁_𝑛=0𝐞̂^𝑛+1−𝐞̂^𝑛2+𝜈Δ𝑡^𝑛=𝑁_𝑛=0𝐞̃^𝑛+12₁+𝐞̂^𝑛+12₁≤𝑀Δ𝑡². (5.46)

Proof. Taking the inner product of (5.27) with 2Δ𝑡𝐞̂^𝑛+1, we have 𝐞̂^𝑛+12−𝐞̂^𝑛2+𝐞̂^𝑛+1−𝐞̂^𝑛2+𝜈Δ𝑡𝐞̃^𝑛+12₁+𝜈Δ𝑡𝐞̂^𝑛+12₁=2Δ𝑡𝑏𝐮^𝑛,𝐮^𝑛+1,𝐞̂^𝑛+1−2Δ𝑡𝑏𝐮𝑡_𝑛+1,𝐮𝑡_𝑛+1,𝐞̂^𝑛+1+2Δ𝑡𝐑^𝑛,𝐞̂^𝑛+1. (5.47) The estimates below are obtained on the right-hand term of (5.47): 2Δ𝑡𝐑^𝑛,𝐞̂^𝑛+1≤𝜈Δ𝑡𝐞̂10^𝑛+121+𝑐Δ𝑡^2𝑡𝑛+1_𝑡𝑛𝐮_𝑡𝑡²₋₁𝑑𝑡. (5.48) For the nonlinear term, similar to (5.6), the below estimates are obtained: −2Δ𝑡𝑏𝑢𝑡_𝑛+1,𝐞̃^𝑛+1,𝐞̂^𝑛+1≤𝜈Δ𝑡𝐞̂10^𝑛+121+𝑐Δ𝑡𝐞̃^𝑛+12. (5.49) For the remainder of nonlinear term, it yields 2Δ𝑡𝑏𝐮𝑡_𝑛−𝐮𝑡_𝑛+1,𝐞̃^𝑛+1,𝐞̂^𝑛+1≤𝐮𝑡_𝑛−𝐮𝑡_𝑛+11𝐞̃^𝑛+1₁𝐞̂^𝑛+1₁≤𝑐Δ𝑡^2𝑡𝑛+1_𝑡𝑛𝐮_𝑡²₁𝑑𝑡+𝜈Δ𝑡𝐞̂10^𝑛+121,−2Δ𝑡𝑏𝐞̂^𝑛,𝐮^𝑛+1,𝐞̂^𝑛+1=2Δ𝑡𝑏𝐞̂^𝑛,𝐞̃^𝑛+1,𝐞̂^𝑛+1−2Δ𝑡𝑏𝐞̂^𝑛,𝐮𝑡_𝑛+1,𝐞̂^𝑛+1. (5.50) Thus, 2Δ𝑡𝑏𝐞̂^𝑛,𝑒̃^𝑛+1,𝐞̂^𝑛+1≤𝑐Δ𝑡𝐞̂^𝑛₁𝐞̂^𝑛+1₁𝐞̃^𝑛+11/2𝐞̃^𝑛+1₁^1/2≤𝑐Δ𝑡𝐞̂^𝑛₁𝐞̂^𝑛+1₁𝐞̃^𝑛+11/2≤𝑐Δ𝑡^5/4𝐞̂^𝑛₁𝐞̂^𝑛+1₁≤𝑐Δ𝑡^3/2𝐞̂^𝑛2₁+𝜈Δ𝑡𝐞̂10^𝑛+121,−2Δ𝑡𝑏𝐞̂^𝑛,𝐮𝑡_𝑛+1,𝐞̂^𝑛+1≤𝑐Δ𝑡𝐞̂^𝑛₀𝐮𝑡_𝑛+12𝐞̂^𝑛+1₁≤𝑐Δ𝑡𝐞̂^𝑛₀𝐞̂^𝑛+1₁≤𝜈Δ𝑡𝐞̂10^𝑛+121+𝑐Δ𝑡𝐞̂^𝑛2, (5.51) where we have used Theorem 5.1 and formula (2.6).
Taking the sum of the formula (5.47) for 𝑛 from 0 to 𝑁, together with the above estimates, we get 𝐞̂^𝑁+12+^𝑛=𝑁_𝑛=0𝐞̂^𝑛+1−𝐞̂^𝑛2+𝜈Δ𝑡^𝑛=𝑁_𝑛=0𝐞̃^𝑛+12₁+12𝑒̂^𝑛+12₁=𝐶Δ𝑡^2𝑇₀𝐮_𝑡𝑡²₋₁𝑑𝑡+𝑐Δ𝑡^𝑛=𝑁_𝑛=0𝐞̃^𝑛+12+𝑐Δ𝑡^2𝑇₀𝐮_𝑡²₁𝑑𝑡+𝑐Δ𝑡^3/2𝑛=𝑁_𝑛=0𝐞̂^𝑛2₁+𝑐Δ𝑡^𝑛=𝑁_𝑛=0𝐞̂^𝑛2. (5.52) By virtue of the formula (5.23) and the regularity assumption (A2), (A3), we obtain𝐞̂^𝑁+12+^𝑛=𝑁_𝑛=0𝐞̂^𝑛+1−𝐞̂^𝑛2+𝜈Δ𝑡^𝑛=𝑁_𝑛=0𝐞̃^𝑛+12₁+12𝐞̂^𝑛+12₁=𝐶Δ𝑡²+𝑐Δ𝑡^3/2𝑛=𝑁_𝑛=0𝐞̂^𝑛2₁+𝑐Δ𝑡^𝑛=𝑁_𝑛=0𝐞̂^𝑛+12. (5.53) For sufficiently small Δ𝑡, we can take the last term to left side and apply the discrete Gronwall lemma to the last inequality, so the proof is completed.

Theorem 5.5. Under the regularity assumptions (A1)–(A3), for small enough Δ𝑡, there exists a positive constant 𝛿, such that 𝛿<1/8 and 4𝑀√Δ𝑡‖𝑒^𝑁̂‖²𝜈⁻¹<1, then 𝑒̃^𝑁_𝑛=0^𝑛+1−𝑒̃^𝑛2≤𝑚Δ𝑡². (5.54)

Proof. We shift the index 𝑛+1→𝑛 in (3.5) to get 𝐞̂^𝑛−𝐞̃^𝑛+1Δ𝑡2𝜈Δ𝐮^𝑛−∇𝑝^𝑛−2𝜔×𝐮^𝑛+𝐟𝑡_𝑛=0, (5.55) and taking the sum with (5.5), we obtain 𝐞̃^𝑛+1−𝐞̃^𝑛−𝜈Δ𝑡2Δ𝐞̃^𝑛+1=𝐮^𝑛⋅∇𝐮^𝑛+1−𝐮𝑡_𝑛+1⋅∇𝐮𝑡_𝑛+1+𝐑^𝑛−∇𝑝𝑡_𝑛+1−𝑝^𝑛+𝐟𝑡_𝑛+1−𝐟𝑡_𝑛+𝜈2Δ𝐞̂^𝑛+𝜈2Δ𝐮𝑡_𝑛+1−𝐮𝑡_𝑛−2𝜔×𝐮𝑡_𝑛+1−𝐮𝑡_𝑛−2𝜔×𝐞̂^𝑛. (5.56) Taking the inner product of (5.56) with Δ𝑡(𝐞̃^𝑛+1−𝐞̃^𝑛), the left-hand term of (5.56) can be written as 𝐞̃^𝑛+1−𝐞̃^𝑛2+𝜈Δ𝑡4𝐞̃^𝑛+12₁+𝐞̃^𝑛+1−𝑒̃^𝑛2₁−𝐞̃^𝑛2₁. (5.57) Now, we give the estimates of the right-hand term of (5.56): −2Δ𝑡𝜔×𝐮𝑡_𝑛+1−𝐮𝑡_𝑛,𝐞̃^𝑛+1−𝐞̃^𝑛≤𝑐Δ𝑡^3𝑡𝑛+1_𝑡𝑛𝐮_𝑡²𝑑𝑡+𝛿𝐞̃^𝑛+1−𝐞̃^𝑛2,−2Δ𝑡𝜔×𝐞̂^𝑛,𝐞̃^𝑛+1−𝐞̃^𝑛≤𝑐Δ𝑡²𝐞̂^𝑛2+𝛿𝐞̃^𝑛+1−𝐞̃^𝑛2≤𝑀Δ𝑡³+𝛿𝐞̃^𝑛+1−𝐞̃^𝑛2, (5.58) where we have used Theorem 5.1. Simultaneously, 𝜈Δ𝑡2Δ𝐮𝑡_𝑛+1−𝐮𝑡_𝑛,𝐞̃^𝑛+1−𝐞̃^𝑛=−𝜈Δ𝑡2∇𝐮𝑡_𝑛+1−𝐮𝑡_𝑛,∇𝐞̃^𝑛+1−𝐞̃^𝑛≤𝑐Δ𝑡𝐮𝑡_𝑛+1−𝐮𝑡_𝑛1𝐞̃^𝑛+1−𝐞̃^𝑛₁≤𝑐Δ𝑡𝐮𝑡_𝑛+1−𝐮𝑡_𝑛²₁+𝜈Δ𝑡8𝐞̃^𝑛+1−𝐞̃