Abstract
We present new error estimates for the finite volume and finite difference methods applied to the compressible Navier–Stokes equations. The main innovative ingredients of the improved error estimates are a refined consistency analysis combined with a continuous version of the relative energy inequality. Consequently, we obtain better convergence rates than those available in the literature so far. Moreover, the error estimates hold in the whole physically relevant range of the adiabatic coefficient.
We’re sorry, something doesn't seem to be working properly.
Please try refreshing the page. If that doesn't work, please contact support so we can address the problem.
Avoid common mistakes on your manuscript.
1 Introduction
The Navier–Stokes equations governing the motion of viscous compressible fluids have numerous applications in engineering, physics, meteorology or biomedicine. In this paper we consider the viscous barotropic fluid endowed, for simplicity, with the isentropic pressure–density state equation \(p= a \varrho ^{\gamma },\) where \(a > 0\) is a positive constant, and \(\gamma >1\) denotes the adiabatic coefficient. The global-in-time existence of weak solutions is known for any \(\gamma > \frac{d}{2}\) in the d-dimensional setting, see Lions [20] and Feireisl [7]. More recently, Plotnikov and Vaigant [27] extended the existence theory for any \(\gamma \ge 1\) if \(d=2\). Unfortunately, the multilevel approach used in the existence proof is rather difficult to adapt directly to a numerical scheme; whence the numerical analysis of the problem remains rather incomplete.
In the last few decades, many efficient and robust numerical methods have been proposed to simulate the motion of viscous compressible fluid flows. We refer the reader to the monographs by Dolejší and Feistauer [2], Eymard et al. [5], Feistauer [3], Feistauer et al. [4], Toro [28], and the references therein. Despite a good agreement of the obtained results with experiments, a rigorous convergence analysis with the associated error estimates has been performed only in a few particular cases.
In his truly pioneering work, Karper [18], see also [9], showed convergence (up to a subsequence) of a mixed finite element-finite volume (or discontinuous Galerkin) approximation to a weak solution of the compressible multidimensional Navier–Stokes system under the technical restriction \(\gamma > 3\). His proofs basically follow step by step the existence theory developed in [7] and as such is difficult to adapt to other numerical methods. Moreover, as the weak solutions are not known to be unique, the result holds up to a subsequence and no convergence rate is available.
Recently, see [10,11,12], we have developed a new approach based on the concept of more general dissipative weak (dissipative measure-valued) solutions, which, combined with the weak–strong uniqueness and conditional regularity results, yields a rigorous proof of convergence for the mixed finite element-finite volume, finite volume and finite difference Marker-and-Cell (MAC) methods for any \(\gamma >1\) as long as the sequence of the numerical solution remains uniformly bounded and/or if the strong solution exists. The aim of the present paper is to derive error estimates for the finite volume and the MAC methods for the full range of the adiabatic coefficient \(\gamma >1.\)
There are several results concerning error estimates for the compressible Navier–Stokes equations. Under the assumption of the \(L^2\)-bounds of the discrete derivatives of the numerical solutions, Jovanović [17] studied the convergence rate of a finite volume-finite difference method to the barotropic Navier–Stokes system. In [21, 22] Liu analyzed the errors for \(P^k\) conforming finite element method, \(k\ge 2\), assuming the existence of a suitably regular smooth solution. However, the stability of the method with respect to the discrete energy was not investigated.
Furthermore, Gallouët et al. [14, 15] analyzed the unconditional convergence rates of the mixed finite volume-finite element method [9] and the MAC scheme for \(\gamma >3/2\) in the dimension \(d=3\). Similar results have been obtained by Mizerová and She [24]. All the above mentioned convergence results are based on a discrete version of the relative energy inequality estimating the error between the numerical and the strong solution. The obtained convergence error is \({\mathcal {O}}(h^A),\) where \(h>0\) is a mesh parameter and \(A = \min \left\{ \frac{2\gamma -3}{\gamma },\frac{1}{2} \right\} \), cf. [14, 15, 24]. In particular, the convergence order tends to zero when \(\gamma \rightarrow \frac{3}{2}\) and remains positive only if \(\gamma >3/2.\) Moreover, if \(\gamma \ge 2 \), the convergence rate is only \(\frac{1}{2}\) in the energy norm, though the numerical experiments indicate the second-order convergence rate.
In view of the existing results, the main novelty of the present paper is twofold:
-
Extending the error analysis to the full range \(\gamma > 1\).
-
Improving the convergence rate via a detailed consistency and error analysis.
Following the strategy proposed in the monograph [12, Chapter 9], we combine the standard consistency errors with the “continuous” form of the relative energy inequality. In contrast with the existing methods based on ad hoc construction of an approximate relative energy inequality, the new approach is rather versatile and free of additional discretization errors. In particular, we can handle any consistent energy stable numerical method in the same fashion. We focus on the finite volume method proposed in [12] and the MAC method from [24]. The application to the mixed finite element-finite volume method of Karper [18] was studied independently and presented in the recent work by Novotný and Kwon [19]. Compared to the previous results of Gallouët et al. [14, 15], we employ the consistency formulation of the numerical solution where the test function is smooth. This new approach avoids the complicated integration by parts formulae on the discrete level and improves the convergence rates of the MAC method presented in [14, 24].
The paper is organized in the following way. After presenting the continuous model and the corresponding relative energy, we formulate the numerical schemes: the finite volume and the MAC method, see Sect. 2. Next, we discuss their energy stability and consistency. The main results on the error estimates are formulated and proved in Sect. 3.
1.1 Compressible Navier–Stokes system
We begin with formulating the compressible Navier–Stokes system
in the time-space cylinder \([0,T] \times \Omega \), \(\Omega \subset R^d, d=2,3\), where \(\varrho \) is the density, \(\varvec{u}\) is the velocity field, and \( \mathbb {S}\) is the viscous stress tensor given by
The pressure is assumed to satisfy the isentropic law
To avoid technical problems related to a proper numerical approximation of the physical boundary, we impose the periodic boundary conditions and identify the computational domain with the flat torus
The system (1.1) is supplemented with finite energy initial data \((\varrho _0, \varvec{u}_0): \mathbb {T}^d\rightarrow \mathbb {R}^+ \times \mathbb {R}^d\),
where \( {P} \) is the so-called pressure potential, \( {P} (\varrho ) = \frac{a \varrho ^{\gamma }}{\gamma -1}\) for the isentropic gas law (1.2).
1.2 Relative energy
The main tool to evaluate the distance between numerical and strong solutions is the relative energy functional, cf. [8]:
As pointed out, relative energy functionals are often used to estimate the distance between a suitable weak solution and the strong solution; whence yielding the weak-strong uniqueness property. Recently, a discrete version of the relative energy has been applied in the error analysis of numerical schemes, see [14, 15, 24].
1.3 Classical solutions
It will be useful to identify the regularity class of smooth (classical) solutions to the Navier–Stokes system (1.1) inherited from the initial data (1.3). The following result can be deduced from [1, Theorem 3.3] and [6, Proposition 2.2].
Proposition 1.1
Let the initial data belong to the class
Let \((\varrho , \varvec{u})\) be a weak solution to problem (1.1) originating from the initial data (1.3) such that
Then \((\varrho , \varvec{u})\) is a classical solution of (1.1)–(1.3) in \([0,T] \times \mathbb {T}^d\).
If, in addition, \(\varrho _0\), \(\varvec{u}_0\) belong to the class
then \(\varrho \in C([0,T]; W^{k,2}(\mathbb {T}^d))\), \(\varvec{u}\in C([0,T]; W^{k,2}(\mathbb {T}^d; \mathbb {R}^d))\), and the following estimate hold
where D depends solely on \(T, {\bar{r}}, {\bar{u}}\) and the initial data \((\varrho _0, \varvec{u}_0)\) via the norm \(\Vert (\varrho _0, \varvec{u}_0) \Vert _{W^{k,2}(\mathbb {T}^d; \mathbb {R}^{d + 1})}\) and \(\min _{x \in {\mathbb {T}^d}} \varrho _0(x).\)
Proof
The first part was proved in [6, Proposition 2.2] via the local existence theory by Valli and Zajaczkowski [26] combined the weak–strong uniqueness principle and the conditional regularity result by Sun et al. [25]. In particular, the bounds (1.6) were established for \(k = 3\), \(\ell = 1\).
Next, as shown in [1, Theorem 3.3], the solution inherits higher Sobolev regularity from the data as long as the norm \(\Vert \varvec{u}\Vert _{C([0,T]; W^{2,\infty }(\mathbb {T}^d; \mathbb {R}^d))}\) is controlled. In particular, the estimates (1.6) can be established. Similarly to Gallagher [13], the proof in [1] is based on the particular isentropic form of the pressure that enables to transform the problem to a parabolic perturbation of a symmetric hyperbolic system. \(\square \)
2 Numerical methods
First, we introduce suitable notation. By c we denote a positive constant independent of the discretization parameters \(\Delta t\) and h. We shall frequently write \(A \lesssim B\) if \(A \le cB\) and \(A\approx B\) if \(A \lesssim B\) and \(B \lesssim A\). We also write \(c\in \text {co}\{a,b\} \) if \(\min (a,b)\le c\le \max (a,b)\). Moreover, we denote by \(\Vert \cdot \Vert _{L^p}\), \(\Vert \cdot \Vert _{L^pL^q}\), and \(\Vert \cdot \Vert _{L^pW^{q,s}}\) the norms \(\Vert \cdot \Vert _{L^p(\mathbb {T}^d)},\) \(\Vert \cdot \Vert _{L^p(0,T;L^q(\mathbb {T}^d))}\), and \(\Vert \cdot \Vert _{L^p(0,T;W^{q,s}(\mathbb {T}^d))},\) respectively.
2.1 Time discretization
We divide the time interval [0, T] into \(N_t\) equidistant parts with a fixed time increment \(\Delta t\) (\( =T/ N_t\)). For a function \(f^n\) given at the discrete time instances \(t_n = n \Delta t\), \(n=0,1,\ldots ,N_t\), we define a piecewise constant approximation f(t) in the following way
The time derivative is approximated by the backward Euler method
2.2 Space discretization
To begin, we introduce a uniform structured mesh including primary, dual, and bidual grids.
Primary grid
We call \({\mathcal {T}}\) the primary grid with the following properties and notations:
-
The domain \(\mathbb {T}^d\) is divided into compact uniform quadrilaterals \( \mathbb {T}^d=\bigcup _{K\in {\mathcal {T}}} K\), where \({\mathcal {T}}\) is the set of all elements that form the primary grid.
-
\({\mathcal {E}}\) denotes the set of all faces of the primary grid \({\mathcal {T}}\). Given an element \(K\in {\mathcal {T}}\), \({\mathcal {E}}(K)\) is the set of its faces; \({\mathcal {E}}_i\) is the set of all faces that are orthogonal to the unit basis vector \(\varvec{e}_i\); \( {\mathcal {E}}_i(K)= {\mathcal {E}}(K)\cap {\mathcal {E}}_i\) for any \(i \in \{1, \ldots , d \}\).
-
h denotes the uniform size of the grid, meaning \(|x_K-x_L|= h\) for any neighbouring elements K and L, where \(x_K\) and \(x_L\) are the centers of K and L, respectively.
-
\(\sigma _{K\text {,} i -}\) and \(\sigma _{K\text {,} i +}\) denote the left and right face of an element K in the \( i{\text {th}}\)-direction, respectively.
-
\({\mathcal {N}}(K)\) denotes the set of all neighbouring elements of \(K \in {\mathcal {T}}.\)
-
\(\sigma = K|L\) denotes the face \(\sigma \) that separates the elements K and L. Moreover, \(\sigma = \overrightarrow{K|L}\) means \(\sigma = K|L\) and \(x_L - x_K = h\varvec{e}_i\) for some \(i \in \{1, \ldots , d \}\).
-
\(\varvec{n}\) denotes the outer normal of a generic face \(\sigma \) and \(\varvec{n}_{\sigma , K}\) denotes the outer normal vector to a face \(\sigma \in {\mathcal {E}}(K).\)
Dual grid
The dual of the primary grid is determined as follows.
-
For any face \(\sigma =K|L\in {\mathcal {E}}_i\), a dual cell is defined as \(D_\sigma =D_{\sigma ,K} \cup D_{\sigma ,L}\), where \(D_{\sigma ,K} = \{x\in K, x_i \in \text {co}\{(x_K)_i, (x_\sigma )_i \} \}\), see Fig. 1a for a two dimensional graphic illustration.
-
\({\mathcal {D}}_i=\left\{ D_\sigma \; | \; \sigma \in {\mathcal {E}}_i\right\} \), \(i \in \{1, \ldots , d \}\), represents the \( i{\text {th}}\) dual grid of \({\mathcal {T}}\). Note that for each fixed \(i \in \{1, \ldots , d \}\) it holds
$$\begin{aligned} \mathbb {T}^d=\bigcup _{\sigma \in {{\mathcal {E}}_i}} D_{\sigma }, \;\; \textrm{int}(D_\sigma )\cap \textrm{int}( D_{\sigma '})=\emptyset \text{ for } \sigma ,\sigma '\in {{\mathcal {E}}_i},\,\sigma \ne \sigma '. \end{aligned}$$ -
\({\widetilde{{\mathcal {E}}}}_i\) is the set of all faces of the \( i{\text {th}}\) dual grid \({\mathcal {D}}_i\) and \( {\widetilde{{\mathcal {E}}}}_{i,j}= \{ \epsilon \in {\widetilde{{\mathcal {E}}}}_i | \epsilon \text{ is } \text{ orthogonal } \text{ to } \varvec{e}_j\}\).
-
A generic face of a dual cell \(D_\sigma \) is denoted as \(\epsilon \in {\widetilde{{\mathcal {E}}}}(D_\sigma )\), where \({\widetilde{{\mathcal {E}}}}(D_\sigma )\) denotes the set of all faces of \(D_\sigma .\)
-
\(\epsilon = D_\sigma | D_{\sigma '}\) denotes a dual face that separates the dual cells \(D_\sigma \) and \(D_{\sigma '}\). Moreover, \(\epsilon = \overrightarrow{ D_\sigma | D_{\sigma '} }\) means \(\epsilon = D_\sigma | D_{\sigma '}\) and \(x_{\sigma '} - x_{\sigma } = h \varvec{e}_i\) for some \(i \in \{ 1,\ldots , d\}\).
-
\( {\mathcal {N}}^\star (\sigma )\) denotes the set of all faces whose associated dual elements are the neighbours of \(D_\sigma ,\) i.e.,
$$\begin{aligned} {\mathcal {N}}^\star (\sigma )=\{\sigma '\ |\ D_{\sigma '} \text{ is } \text{ a } \text{ neighbour } \text{ of } D_\sigma \}. \end{aligned}$$
Bidual grid
-
Similarly to the definition of the dual cell, a bidual cell \(D_\epsilon := D_{\epsilon ,\sigma } \cap D_{\epsilon ,\sigma '}\) associated to \(\epsilon = D_\sigma | D_{\sigma '} \in {\widetilde{{\mathcal {E}}}}_{i,j}\) is defined as the union of adjacent halves of \(D_\sigma \) and \(D_{\sigma '}\), where \(D_{\epsilon ,\sigma } = \{x \in D_{\sigma }| x_j \in \text {co}\{(x_\sigma )_j, (x_\epsilon )_j \} \} \) see Fig. 1b for a two dimensional graphic illustration.
-
\({\mathcal {B}}_{i,j}\) denotes the \(j^{\textrm{th}}\) dual grid of \({\mathcal {D}}_i\), which is a set of all bidual cells associated with the bidual faces of \( {\widetilde{{\mathcal {E}}}}_{i,j}\). Note that \({\mathcal {B}}_{i,j}= {\mathcal {T}}\) in the case of \(i=j\).
Discrete function spaces. We introduce the following spaces of piecewise constant functions:
The corresponding projections read
where \(1_K\) and \(1_{D_\sigma }\) are the characteristic functions. Further, for any \(\varvec{\phi }=(\phi _1, \ldots , \phi _d)\) we denote \( \Pi _ {\mathcal {E}}\varvec{\phi }= \left( \Pi _ {\mathcal {E}}^{(1)}\phi _1, \ldots , \Pi _ {\mathcal {E}}^{(d)}\phi _d\right) .\) Moreover, for any bidual grid \(D_\epsilon \) we define
2.3 Discrete operators
Average and jump. First, for a piecewise smooth function \(f_h\), we define its trace
Then for any \(r_h \in Q_h\) we define the average operator
If in addition, \(\sigma \in {\mathcal {E}}_i\) for an \(i\in \{1, \ldots , d\}\), we write as and denote
Analogously to the average operator, we define the jump operator for \(r_h \in Q_h\) as
Further, for vector-valued functions and we define
Note that for any we have .
Gradient operator. For any \(r_h \in Q_h\) and we introduce the following gradient operators.
where
Furthermore, for any and \(\phi \in W^{1,2}(\mathbb {T}^d)\) we set
Here, \( \partial _{\mathcal {T}}^{(i)} \) is defined for any \( u_{i,h} \in W_{i,h}\), \(i \in \{1, \ldots , d \}\) as
Note that for any \(r_h \in Q_h\) and \( u_{i,h} \in W_{i,h}\), there hold
Divergence operator. For \(\varvec{u}_h\in \textbf{W}_h\) and \(\varvec{v}_h\in \textbf{Q}_h\) we define the following discrete divergence operators adjoint to the above discrete gradient operators
It is easy to observe for any that
Upwind flux. Given a velocity field , the upwind flux function for \(r_h \in Q_h\) is given by
where
To approximate nonlinear convective terms we apply the following diffusive upwind flux
For we define a vector-valued upwind flux componentwise
2.4 Preliminary estimates and inequalities
In this section we present preliminary material. First, it is easy to check that the following integration by parts formulae hold, see e.g. [16, Lemma 2.1].
Lemma 2.1
Let \(r_h, \phi _h \in Q_h,\) and \(\varvec{u}_h, \varvec{\phi }_h \in \textbf{W}_h.\) Then
Next, we report the following useful lemmas whose proofs are presented in Appendix A.
Lemma 2.2
For any \(r_h \in Q_h\), \(\varvec{v}_h\in \textbf{Q}_h\), \(\varvec{u}_h\in \textbf{W}_h\), \(\psi \in W^{1,2}(\mathbb {T}^d)\) and \(\varvec{U}\in W^{1,2}(\mathbb {T}^d;\mathbb {R}^d)\), there hold
Lemma 2.3
For any \(\varvec{u}_h\in \textbf{W}_h\), \(\varvec{v}_h\in \textbf{Q}_h\) and \(\psi \in W^{1,2}(\mathbb {T}^d)\) there hold
Lemma 2.4
For any \(\varvec{u}_h\in \textbf{W}_h\), \(\varvec{v}_h\in \textbf{Q}_h\) and \(\varvec{U}\in W^{2,2}(\mathbb {T}^d; \mathbb {R}^d)\), we have
Lemma 2.5
Let \(\varvec{v}_h\in \textbf{Q}_h\), \(\varvec{u}_h\in \textbf{W}_h\), \(\varvec{U}\in W^{2,2}(\mathbb {T}^d; \mathbb {R}^d)\), and \(\mathbf {\Phi }\in W^{3,2}(\mathbb {T}^d; \mathbb {R}^d)\). Then for any \(i,j\in \{1, \ldots , d \}\), we have
2.5 Finite volume and finite difference methods
We proceed by presenting a finite volume and a finite difference numerical method that will be used to approximate the Navier–Stokes system (1.1)–(1.3). Both methods have been already successfully applied in numerical simulations, see, e.g., [12]. In our recent work [11, 12, 24], the convergence was shown for \(\gamma > 1\) via the concept of dissipative measure-valued solutions. However, the error analysis was missing for the finite volume method and suboptimal for the finite difference method.
2.5.1 Finite volume method
We introduce the finite volume (FV) method approximating the Navier–Stokes system (1.1)–(1.3).
Definition 2.6
(FV scheme) Given the initial data (1.3), we set \((\varrho _h^0,\varrho _h^0 \varvec{u}_h^0) =(\Pi _Q\varrho _0, \Pi _Q[\varrho _0 \varvec{u}_0])\). The FV approximation \((\varrho _h^n, \varvec{u}_h^n) \in Q_h\times \textbf{Q}_h,\) \(n=1,\dots , N,\) of the Navier–Stokes system (1.1)–(1.3) is a solution of the following system of algebraic equations:
where \(\nu = \frac{d-2}{d} \mu + \lambda .\)
2.5.2 Finite difference MAC method
We proceed by presenting the finite difference MAC scheme that is based on a staggered grid approach. On the one hand, the discrete density \(\varrho _h\) and pressure \(p_h=p(\varrho _h)\) are approximated on the primary grid \({\mathcal {T}}\). On the other hand, the \( i{\text {th}}\) component of the velocity field \( u_{i,h} \) is approximated on the \( i{\text {th}}\) dual grid \({\mathcal {D}}_i\). The MAC scheme reads as follows.
Definition 2.7
(MAC scheme) Given the initial data (1.3), we consider \((\varrho _h^0,\varrho _h^0 \Pi _Q\varvec{u}_h^0) =(\Pi _Q\varrho _0, \Pi _Q[\varrho _0 \varvec{u}_0]).\) The MAC approximation of the Navier–Stokes system (1.1)–(1.3) is a sequence \( (\varrho _h^n,\varvec{u}_h^n ) \in Q_h\times \textbf{W}_h\), \(n=1,2,\dots ,N\), which solves the following system of algebraic equations:
where \(\nu = \frac{d-2}{d} \mu + \lambda .\)
In what follows, we will denote by \(\varrho _h(t), \varvec{u}_h(t)\) the piecewise constant approximations of \(\varrho _h^n, \varvec{u}_h^n\), \(n=0,1, \dots , N\) on the time interval [0, T], see Sect. 2.1. We note that both methods, the FV method (2.10) as well as the MAC method (2.11), preserve the positivity of density and conserve the mass
where \( M:= \int _{\mathbb {T}^d} \varrho _0 \,\textrm{d}x \) denotes the fluid mass, see e.g. [12, Lemma 11.2].
2.6 Energy stability
The essential feature of any numerical scheme is its stability. We now recall the energy stability of both numerical methods introduced above, see [12, Theorem 11.1 and 14.1]
Lemma 2.8
(Energy estimates) Let \((\varrho _h, \varvec{u}_h)\) be a numerical solution obtained either by the FV scheme (2.10) or by the MAC scheme (2.11) with \(\gamma >1\). Then for all \(\tau \in (0,T)\), it holds
where \( E_0=\int _{\mathbb {T}^d} \Big (\frac{1}{2} \varrho _0 |\varvec{u}_0|^2 + {P} (\varrho _0) \Big ) \,\textrm{d}x\) is the initial energy and
Moreover, there exists \(c>0\) which may depend on the fluid mass M and the initial energy \(E_0\) but is independent of the parameters h and \(\Delta t\) such that
2.7 Consistency formulation
The next important ingredient of our approach is the consistency formulation of the numerical scheme.
Lemma 2.9
(Consistency formulation) Let \((\varrho _h, \varvec{u}_h)\) be either a solution of the FV scheme (2.10) or the MAC scheme (2.11) with \(\Delta t \approx h \in (0,1)\), \(\gamma >1\) and \(\varepsilon >-1\).
Then for all \(\tau \in (0,T)\), \(\phi \in L^\infty (0,T; W^{2,\infty }(\mathbb {T}^d)),\) \( \partial _t ^2 \phi \in L^\infty ((0,T) \times \mathbb {T}^d)\) and \(\varvec{\phi }\in L^\infty (0,T; W^{2,\infty }(\mathbb {T}^d; \mathbb {R}^d))\), \( \partial _t ^2 \varvec{\phi }\in L^\infty ((0,T) \times \mathbb {T}^d; \mathbb {R}^d)\) there hold
where the consistency errors are bounded as follows:
Here, the constant \(C_\varrho \) depends on
and \(C_{\varvec{m}}\) depends on
Further, the exponents \(\beta _D\) and \(\beta _M\) are given by
Remark 1
Consistency formulation for the FV and MAC method was introduced in [12, Theorem 11.2] and [12, Theorem 14.2], respectively. Instead of an abstract consistency error identified in [12], Lemma 2.9 provides an explicit bound in terms of the numerical step and regularity of the associated test function. Moreover, we improve the result of [12] by requiring less regularity of the test functions.
Remark 2
Following the recent result of Lukáčová-Medvid’ová et al. [23, Lemma B.2], we could improve the estimates on \(\beta _D\) and \(\beta _M\).
Proof of Lemma 2.9
The consistency errors arising from the time derivative term can be evaluated in the following way. First, by a direct calculation, we obtain
for any \(\tau \in [t_n, t_{n+1})\), \(n=1,\dots , N_T\), where
Collecting the above estimates we obtain from (2.16) that
whenever \(\tau \in [t_n, t_{n+1})\), where (\(r_h, \varphi \)) stands for \((\varrho _h,\phi )\) or \((\varrho _h\varvec{u}_h, \varvec{\phi })\).
Analogously as in the proofs of [12, Theorem 11.2] and [12, Theorem 14.2], we obtain
where \(e_{\varvec{m}}^*\) is controlled by
In order to derive (2.15a) it suffices to realize that the time integral from \(\tau \) to \(t^{n+1}\) at the right-hand side of (2.18a) is of order \({\mathcal {O}}(\Delta t)\). Indeed
Combining (2.18a) and (2.19) yields (2.15a).
Similarly, to get (2.15b) we need the following estimate
Substituting the above estimate into (2.18b) we obtain (2.15b), which completes the proof. \(\square \)
Lemma 2.10
(Consistency formulation for a bounded numerical solution) Let the assumptions of Lemma 2.9 hold. Moreover, let \(\varrho _h\) and \(\varvec{u}_h\) be uniformly bounded, i.e., there exist positive constants \(\overline{\varrho }\) and \(\overline{u}\) such that
Then for all \(\tau \in (0,T)\), \(\phi \in L^\infty (0,T; W^{2,\infty }(\mathbb {T}^d))\), \( \partial _t ^2 \phi \in L^\infty ((0,T) \times \mathbb {T}^d)\) and \(\varvec{\phi }\in L^\infty (0,T; W^{2,\infty }(\mathbb {T}^d; \mathbb {R}^d))\), \( \partial _t ^2 \varvec{\phi }\in L^\infty ((0,T) \times \mathbb {T}^d; \mathbb {R}^d)\), there holds
where the consistency errors can be bounded as follows
Here, the constant \(C_\varrho \) depends on
and \(C_{\varvec{m}}\) depends on
Proof
We will present the proof for the FV method, the proof for the MAC method is analogous. First, we denote the errors of the inviscid fluxes as
Analogously as in the proof of [12, Theorem 11.3] we get
In view of assumption (2.20) the errors \(e_1\) and \(e_2\) are controlled by
Now, summing up (2.22) and (2.17) with \(r_h =\varrho _h\), and recalling the estimates (2.19) and (2.24) implies (2.21a). Moreover, summing up (2.23) and (2.17) with \(r_h =\varrho _h\varvec{u}_h\) we get
where \(e_2\) is given in (2.23) and \(e_3,e_4,e_5\) read
Recalling (2.8c), (2.9b), (2.9d) and using Hölder’s inequality, we get
Similarly, recalling (2.8d), (2.9a), (2.9c) and Hölder’s inequality, we get
We are left with the estimate of \(e_5\), which can be controlled due to Hölder’s inequality.
Consequently, collecting the estimates of \(e_2\), \(e_3\), \(e_4\) and \(e_5\) we observe that (2.21b) follows from (2.25), which completes the proof. \(\square \)
3 Error estimates
This section is the heart of the paper. We prove the main result – the convergence rates for the FV (2.10) and MAC (2.11) schemes. If, in addition, the numerical solutions are uniformly bounded, the convergence rates can be improved to the first order.
Theorem 3.1
(Convergence rates) Let \(\gamma > 1\) and the initial data \((\varrho _0, \varvec{u}_0)\) satisfy
Suppose that the Navier–Stokes system (1.1) admits a classical solution \((\varrho ,\varvec{u})\) defined on \([0,T] \times \mathbb {T}^d\), with the initial data \((\varrho _0, \varvec{u}_0).\) Further, let \((\varrho _h, \varvec{u}_h)\) be a numerical solution obtained either by the FV scheme (2.10) or by the MAC scheme (2.11) emanating from the projected initial data \((\varrho _h^0, \varvec{u}_h^0)\).
Then there exists a positive number
such that
and
The convergence rate A reads
Here the constants \(\beta _D\) and \(\beta _M\) are given in (2.15d).
Remark 3
Let us discuss the obtained convergence rate \({\mathcal {O}}(h^A)\) for the choice \(\Delta t=h\) and different values of \(\gamma >1\), \(d=2,3.\)
-
For the case \(d=2\), we obtain the following convergence rate A:
-
Let \(\gamma \ge 2\). Then for any \(\varepsilon \ge 0\) both numerical methods have the first-order convergence rate, i.e. \(A=1.\)
-
Let \(\gamma \in (1,2)\). The convergence rates are different for the FV and MAC schemes.
-
* \(A_{FV}= \min \left\{ 1- \frac{5+3\varepsilon }{6\gamma }, 1, 1+\varepsilon \right\} \). Choosing the optimal value of \(\varepsilon \), \(\varepsilon = - \frac{5}{3+6\gamma }\in (-\frac{5}{9}, -\frac{1}{3})\), the convergence rate \(A_{FV}=1 + \varepsilon \) varies between \(\frac{4}{9}\) for \(\gamma \searrow 1\) and \(\frac{2}{3}\) for \(\gamma \nearrow 2.\)
-
* \(A_{MAC}= \min \left\{ 1- \frac{5+3\varepsilon }{6\gamma }, 1, 1+\varepsilon , 1+\varepsilon -\frac{5+3\varepsilon }{6\gamma } \right\} \) reaches its maximum value \(\frac{6\gamma -5}{6\gamma }>0\) at \(\varepsilon =0\). Thus, the convergence rate varies between \(\frac{1}{6}\) for \(\gamma \searrow 1\) and \(\frac{7}{12}\) for \(\gamma \nearrow 2.\)
-
-
-
For the case \(d=3\), we obtain the following convergence rate A:
-
Let \(\gamma \ge 3\). Then for any \(\varepsilon \ge 0\) both methods have first-order convergence rates, i.e. \(A=1\).
-
Let \(\gamma \in [2,3)\). Then for any \(\varepsilon \ge \frac{2\gamma -3}{\gamma }\) we have \(A= \frac{2\gamma -3}{\gamma }\) and the convergence rate varies between \(\frac{1}{2}\) for \(\gamma = 2\) and 1 for \(\gamma \nearrow 3\).
-
Let \(\gamma \in (1,2)\).
-
* \(A_{FV}= \min \left\{ 1- \frac{2+\varepsilon }{2\gamma }, 1, 1+\varepsilon \right\} \). Choosing an optimal value of \(\varepsilon ,\) \(\varepsilon = - \frac{ 2}{1+2\gamma }\in (-\frac{2}{3}, -\frac{2}{5})\), \(A_{FV} = 1+\varepsilon \) and varies between \(\frac{1}{3}\) for \(\gamma \searrow 1\) and \(\frac{3}{5}\) for \(\gamma \nearrow 2.\)
-
* \(A_{MAC}= \min \left\{ 1- \frac{2+\varepsilon }{2\gamma }, 1, 1+\varepsilon , 1+\varepsilon -\frac{2+\varepsilon }{2\gamma } \right\} \) reaches its maximum value \(\frac{\gamma -1}{\gamma } > 0\) at \(\varepsilon =0\). Note that \(A_{MAC}\) varies between 0 when \(\gamma \searrow 1\) and \(\frac{1}{2}\) when \(\gamma \nearrow 2.\)
-
-
Remark 4
In view of the above results, the convergence rates available in the literature, see e.g. [14, 15, 24], are not optimal. Indeed, for \(d=3\) and \(\gamma =\frac{3}{2}\), they degenerate to 0. Moreover, no error analysis is available for \(\gamma < \frac{3}{2}.\) Our approach yields error estimates also for \(\gamma \in (1, \frac{3}{2}].\) In addition, we have better convergence rates, e.g., for \(d=3\) and \(\gamma =\frac{3}{2}\), where the convergence errors are \({\mathcal {O}}(h^{\frac{3}{4} })\) and \({\mathcal {O}}(h^{\frac{1}{3} })\) for the FV and MAC schemes, respectively.
Proof of Theorem 3.1
First, by a straightforward but lengthy calculation, see Appendix D, we observe the following relative energy inequality
Next, we observe the following identities
Then by substituting the above equalities into (3.5) and denoting
we obtain
As \((\varrho ,\varvec{u})\) satisfies the Navier–Stokes system (1.1), we know that
Substituting this equality into (3.6) we get
Rearranging the terms on the right-hand side, we arrive at
where the integrals \(R^E_i, i=1,\ldots ,5\), read
Next, for \(i=1,\ldots , 5\) we analyze \(R^E_i\) such that it can be controlled either by the relative energy or the mesh parameter h.
Term \(R^E_1\). Applying Hölder’s inequality and Lemma B.4 we obtain
where \(C_0^*>0\) depends on \(\Vert \varvec{u}\Vert _{L^\infty W^{2,\infty }}\), \(\Vert \varrho \Vert _{C([0,T] \times \mathbb {T}^d)}, M, E_0, \gamma \), \(\delta \), and \({\underline{r}} = \min _{[0,T] \times \mathbb {T}^d} \varrho \);
\( C_1^*>0\) depends on \(\Vert \varvec{u}\Vert _{L^\infty W^{2,\infty }}\), \( M, E_0\), and \(\gamma \);
\(C_2^*>0\) depends on \(\Vert \varvec{u}\Vert _{L^\infty W^{2,\infty }}\), M, \(E_0, \gamma \), and \(\Vert \nabla _x\varvec{u}\Vert _{L^\infty ( (0,T)\times \mathbb {T}^d)} \).
Term \(R^E_2\). Thanks to Hölder’s inequality we observe the following estimate.
where C depends on \(\Vert \nabla _x\varvec{u}\Vert _{L^\infty ( (0,T)\times \mathbb {T}^d)} \).
Term \(R^E_3\). We analyze the third term \(R^E_3\) in two cases.
First, we consider the case of the FV scheme. In this case \(\varvec{u}_h\in \textbf{Q}_h\), \( \nabla _h\varvec{u}_h= \nabla _{{{\mathcal {D}}}}\varvec{u}_h\) and \(\Pi _Q\varvec{u}_h= \varvec{u}_h\). Thus,
where we have used the equality (2.4), the integration by parts formula (2.3a), and the estimate (2.9b).
Second, we consider the case of \(\varvec{u}_h\in \textbf{W}_h\) obtained by the MAC scheme. In this case, \( \nabla _h\varvec{u}_h= \nabla _{{{\mathcal {B}}}}\varvec{u}_h\) and the term \(R^E_3\) can be estimated in the following way
where we have applied (2.8a), Hölder’s inequality and the estimate (2.9b).
Consequently, we have for both cases
where the constant C depends on \(\mu ,\) the initial energy \(E_0\) and \(\Vert \varvec{u}\Vert _{L^2 W^{2,2}}\).
Term \(R^E_4\). We analyze the term \(R^E_4\) also in two cases.
First, for \(\varvec{u}_h\in \textbf{Q}_h\) obtained by the FV method, we have \( \textrm{div}_h\varvec{u}_h= \textrm{div} _{\mathcal {T}}^{Q}\varvec{u}_h\), \(\Pi _Q\varvec{u}_h=\varvec{u}_h\) and thus
where we have used the identity (2.8d), Hölder’s inequality, the estimates (2.9a) and (2.9c).
Second, for the case of \(\varvec{u}_h\in \textbf{W}_h\) we have
where (2.8b), Hölder’s inequality, the estimates (2.9a) and (2.9c) were applied.
Consequently, we have for both cases
where the constant C depends on \(\nu \), initial energy \(E_0\), and \(\Vert \varvec{u}\Vert _{L^2 W^{2,2}}\).
Term \(R^E_5\). The estimate of \(R^E_5\) is straightforward by applying Hölder’s inequality, i.e.,
where C depends on \(\Vert \textrm{div} _x\varvec{u}\Vert _{L^\infty ((0,T)\times \mathbb {T}^d)}\).
Consequently, collecting the above estimates of \(R^E_i\) for \(i=1,\cdots , 5\), we find
Applying the standard projection error estimates we get
where C depends on \(\Vert \varrho _0\Vert _{C}\) and \(\Vert \varvec{u}_0\Vert _{L^2 W^{1,2}}\).
Consequently, by choosing \(\delta < \frac{\mu }{C_1^*}\), substituting (3.8) into (3.7), using Gronwall’s lemma and recalling the consistency error (2.15c), we may infer that
for \(\Delta t< \frac{1}{C_0^*}\). Here, the constant C depends on \(\Vert \varrho \Vert _{L^{\infty } W^{2,\infty } }, \Vert \varvec{u}\Vert _{L^{\infty } W^{2,\infty } } \) and the exponent A is given by (3.4).
Finally, we combine the above estimate with Lemmas C.1 and B.2 in order to obtain (3.2) and (3.3), respectively. Note that \(E_0\) and M are bounded by the norm \(\Vert (\varrho _0, \varvec{u}_0)\Vert _{W^{k,2}(\mathbb {T}^d; R^{d + 1})}.\) Due to Proposition 1.1 all terms depending on the norms of the exact solution \((\varrho ,\varvec{u})\) as well as \({\underline{r}}\) are bounded by a constant \(c=c( T, \Vert (\varrho _0, \varvec{u}_0)\Vert _{W^{k,2}(\mathbb {T}^d; R^{d + 1})}, \Vert (\varrho , \varvec{u})\Vert _{C([0,T]\times \mathbb {T}^d; R^{d + 1})})\) which finishes the proof. \(\square \)
Finally, we observe that under the assumption that the numerical solutions \((\varrho _h,\varvec{u}_h)\) are uniformly bounded, the above error estimates can be improved. Indeed, applying Lemmas 2.10, C.1, and B.2 we derive the first-order error rate.
Theorem 3.2
(Error rates for bounded numerical solutions) In addition to the hypotheses of Theorem 3.1, let the numerical solution \((\varrho _h, \varvec{u}_h)\) be uniformly bounded,
Then there exists a positive number
such that
for all \(\tau \in [0,T]\), and
4 Conclusion
In this paper we have presented improved error estimates for two well-known numerical methods applied to compressible Navier–Stokes equations. Specifically, we consider the upwind finite volume method and the Marker-and-Cell (MAC) method with implicit time discretization and piecewise constant approximation in space. However, the approach presented in the paper can be applied also to other well-known numerical methods for compressible Navier–Stokes equations.
The novelty of our approach lies in the use of a continuous form of relative energy inequality combined with a refined consistency analysis. Thus, following the framework of the Lax equivalence theorem it suffices to show the (energy) stability, cf. Lemma 2.8, and the consistency of a numerical scheme, cf. Lemma 2.9, in order to obtain the convergence rates for the scheme. Indeed, consistency errors directly yield global errors in the relative energy. To obtain the corresponding error estimates we only assume that the initial data are sufficiently regular and that a strong solution exists globally in time. The error estimates presented in Theorem 3.1 improve the results already presented in the literature [14, 15, 24], see Remarks 2, 3 for a detailed discussion. In particular, our error estimates hold for the full range of the adiabatic coefficient \(\gamma > 1.\)
Moreover, we have considered a natural hypothesis on uniformly bounded numerical solutions and proved that the error estimates can be further improved, cf. Theorem 3.2. Indeed, we prove that both numerical methods converge with the first-order in time and mesh parameter in terms of the relative energy and with the half order in the \(L^\infty (0,T; L^2(\mathbb {T}^d))\)-norm for the density and momentum, as well as in the \(L^2((0,T) \times \mathbb {T}^d)\)-norm for the velocity.
References
Breit, D., Feireisl, E., Hofmanová, M.: Local strong solutions to the stochastic compressible Navier–Stokes system. Commun. Partial Differ. Equ. 43(2), 313–345 (2018)
Dolejší, V., Feistauer, M.: Discontinuous Galerkin Method. Springer Series in Computational Mathematics, vol. 48. Springer (2015)
Feistauer, M.: Mathematical Methods in Fluid Dynamics. Pitman Monographs and Surveys in Pure and Applied Mathematics, vol. 67. Longman Scientific & Technical, Harlow (1993)
Feistauer, M., Felcman, J., Straškraba, I.: Mathematical and Computational Methods for Compressible Flow. The Clarendon Press, Oxford University Press (2003)
Eymard, R., Gallouët, T., Herbin, R.: Finite volume methods. Handb. Numer. Anal. 7, 713–1018 (2000)
Feireisl, E., Hošek, R., Maltese, D., Novotný, A.: Error estimates for a numerical method for the compressible Navier–Stokes system on sufficiently smooth domains. ESAIM Math. Model. Numer. Anal. 51(1), 279–319 (2017)
Feireisl, E.: Dynamics of Viscous Compressible Fluids. Oxford Lecture Series in Mathematics and its Applications, Oxford University Press (2004)
Feireisl, E., Jin, B.J., Novotný, A.: Relative entropies, suitable weak solutions, and weak strong uniqueness for the compressible Navier–Stokes system. J. Math. Fluid Mech. 14(4), 717–730 (2012)
Feireisl, E., Karper, T.G., Pokorný, M.: Mathematical Theory of Compressible Viscous Fluids: Analysis and Numerics. Birkhäuser-Verlag, Basel (2017)
Feireisl, E., Lukáčová-Medvid’ová, M.: Convergence of a mixed finite element-finite volume scheme for the isentropic Navier–Stokes system via the dissipative measure-valued solutions. Found. Comput. Math. 18(3), 703–730 (2018)
Feireisl, E., Lukáčová-Medvid’ová, M., Mizerová, H., She, B.: Convergence of a finite volume scheme for the compressible Navier–Stokes system. ESAIM: M2AN 53(6), 1957–1979 (2019)
Feireisl, E., Lukáčová-Medvid’ová, M., Mizerová, H., She, B.: Numerical Analysis of Compressible Fluid Flows. Springer (2021)
Gallagher, I.: A remark on smooth solutions of the weakly compressible periodic Navier–Stokes equations. J. Math. Kyoto Univ. 40(3), 525–540 (2000)
Gallouët, T., Maltese, D., Novotný, A.: Error estimates for the implicit MAC scheme for the compressible Navier–Stokes equations. Numer. Math. 141, 495–567 (2019)
Gallouët, T., Herbin, R., Maltese, D., Novotný, A.: Error estimates for a numerical approximation to the compressible barotropic Navier–Stokes equations. IMA J. Numer. Anal. 36(2), 543–592 (2016)
Hošek, R., She, B.: Stability and consistency of a finite difference scheme for compressible viscous isentropic flow in multi-dimension. J. Numer. Math. 26(3), 111–140 (2018)
Jovanović, V.: An error estimate for a numerical scheme for the compressible Navier–Stokes system. Kragujev. J. Math. 30, 263–275 (2007)
Karper, T.: A convergent FEM-DG method for the compressible Navier–Stokes equations. Numer. Math. 125(3), 441–510 (2013)
Kwon, Y., Novotný, A.: Consistency, convergence and error estimates for a mixed finite element-finite volume scheme to compressible Navier–Stokes equations with general inflow/outflow boundary data. IMA J. Numer. Anal. 42(1), 107–164 (2022)
Lions, P.L.: Mathematical Topics in Fluid Mechanics. Vol. 2: Compressible Models. Oxford University Press (1998)
Liu, B.: The analysis of a finite element method with streamline diffusion for the compressible Navier–Stokes equations. SIAM J. Numer. Anal. 38, 1–16 (2000)
Liu, B.: On a finite element method for three-dimensional unsteady compressible viscous flows. Numer. Methods Partial Differ. Equ. 20, 432–449 (2004)
Lukáčová-Medvid’ová, M., She, B., Yuan, Y.: Convergence and error estimates of a penalization finite volume method for the compressible Navier–Stokes system. arXiv preprint arXiv:2209.02344
Mizerová, H., She, B.: Convergence and error estimates for a finite difference scheme for the multi-dimensional compressible Navier–Stokes system. J. Sci. Comput. 84(1), 25 (2020)
Sun, Y., Wang, C., Zhang, Z.: A Beale–Kato–Majda blow-up criterion for the 3-D compressible Navier–Stokes equations. J. Math. Pures. Appl. 95(1), 36–47 (2011)
Valli, A., Zajaczkowski, M.: Navier–Stokes equations for compressible fluids: global existence and qualitative properties of the solutions in the general case. Commun. Math. Phys. 103, 259–296 (1986)
Plotnikov, P.I., Weigant, W.: Isothermal Navier–Stokes equations and Radon transform. SIAM J. Math. Anal. 47(1), 626–653 (2015)
Toro, E.F.: Riemann Solvers and Numerical Methods for Fluid Dynamics. A Practical Introduction. Springer (2009)
Funding
Open access publishing supported by the National Technical Library in Prague.
Author information
Authors and Affiliations
Corresponding author
Additional information
Publisher's Note
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
This research was initiated during our “Research in Pairs” stay at the Mathematisches Forschungsinstitut Oberwolfach in 2021.
The research of E.F. and B.S. leading to these results has received funding from the Czech Sciences Foundation (GAČR), Grant Agreement 21-02411S. The Institute of Mathematics of the Academy of Sciences of the Czech Republic is supported by RVO:67985840.
M.L. has been funded by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation)—Project number 233630050—TRR 146 as well as by TRR 165 Waves to Weather. She is grateful to the Gutenberg Research College and Mainz Institute of Multiscale Modelling for supporting her research.
Appendices
Appendix
Proof of the preliminary lemmas
In this section we present the proofs of Lemmas 2.2–2.5.
Proof of Lemma 2.2
First, we calculate
Analogously, we find
which completes the proof.
Proof of Lemma 2.3
First, we calculate
where \(\epsilon ^-\) and \(\epsilon ^+\) are the left and right edges of \(D_\sigma \) in the \( i{\text {th}}\)-direction of the canonical system for \(\sigma \in {\mathcal {E}}_i\). Note that \(D_{\epsilon ^\pm } \subset {\mathcal {T}}\) are elements of the primary grid \({\mathcal {T}}\). Then we can rewrite the above relation as
where we have used (2.3a). This proves (2.6).
The proof of (2.7) follows from (2.5) and (2.3a), specifically,
Proof of Lemma 2.4
First, we recall (2.4) and (2.3a) to derive the first equality
Next, it is easy to check (2.8b) by setting \(\psi = \textrm{div} _x\varvec{U}\) in (2.6), i.e.,
Further, thanks to (2.4) and (2.3a), we observe (2.8c), i.e.,
Finally, by setting \((\left\{ \!\!\left\{ \varvec{v}_h\right\} \!\!\right\} , \textrm{div} _x\varvec{U})\) as \((\varvec{u}_h, \psi )\) into (2.6) we get (2.8d), i.e.
where we have used the identity (2.2).
Proof of Lemma 2.5
Note that the estimates stated in (2.9b)–(2.9d) hold due to the standard interpolation error; whence we omit the proof. Now we prove (2.9a). First, by a direct calculation, we have
where we have used the fact that \( \eth _{{\mathcal {B}}_{i,i}} = \partial _{\mathcal {T}}^{(i)} \) in the last inequality, which proves the first estimate of (2.9a). Analogously, we compute
which proves the second estimate of (2.9a). This concludes the proof of Lemma 2.5.
Sobolev–Poincaré type inequality
First, we recall [12, Theorem 17] for a generalized Sobolev-Poincaré inequality.
Lemma B.1
([12]) For a structure mesh let \(\gamma >1\) and \(\varrho _h\ge 0\) satisfy
where \(\gamma >1\), \(c_M\) and \(c_E\) are positive constants. Then there exists \(c=c(c_M, c_E, \gamma )\) independent of h such that
Now we are ready to show the following lemma.
Lemma B.2
Under the assumption of Lemma B.1 let \((\varrho _h, \varvec{u}_h)\) be a solution obtained either by the FV method (2.10) or the MAC method (2.11). Let \(\varvec{U}\in W^{2,\infty }(\mathbb {T}^d;\mathbb {R}^d)\), then there exists \(C_1 =C_1( M, E_0, \gamma )>0\) and \(C_2=C_2( M, E_0, \gamma , \Vert \nabla _x\varvec{U}\Vert _{L^\infty }, \Vert \varvec{U}\Vert _{W^{2,\infty }} )>0\) such that
where M and \(E_0\) are the fluid mass and initial energy.
Proof
Firstly, by setting \(f_h =\varvec{u}_h-\varvec{U}_h\) for some \(\varvec{U}_h\) belonging to the same discrete space as \(\varvec{u}_h\) in Lemma B.1 we know that
where the constant \(C_1\) depends on \(c_M \equiv M \), \(c_E \equiv E_0\) and \(\gamma \). Note that the choices of \(c_M\) and \(c_E\) are owing to the mass conservation (2.12) and energy stability (2.13).
Next, for \(\varvec{u}_h\in Q_h\) and \(\varvec{u}_h\in \textbf{W}_h\) we set \(\varvec{U}_h= \Pi _Q\varvec{U}\in Q_h\) and \(\varvec{U}_h= \Pi _ {\mathcal {E}}\varvec{U}\in \textbf{W}_h\), respectively. Then by the triangular inequality and projection error we derive
where \(C_2\) depends on \(C_1, \Vert \varvec{U}\Vert _{W^{2,\infty }},\Vert \nabla _x\varvec{U}\Vert _{L^\infty }, M\), and \( \Vert \nabla _x\varvec{U}\Vert _{L^2}\), which proves (B.1).
Finally, we proceed with the proof of (B.2). On the one hand, for the case of \(\varvec{u}_h\in Q_h\) we have \(\Pi _Q\varvec{u}_h=\varvec{u}_h\), meaning (B.2) automatically holds as it is the same as (B.1). On the other hand, for the case of \(\varvec{u}_h\in \textbf{W}_h\) we employ (B.1) and the triangular inequality to derive
where we have used the fact that \(\Vert \textrm{div}_h\varvec{u}_h\Vert _{L^2}^2 \lesssim E_0\) in view of (2.14b), which completes the proof. \(\square \)
Next, we recall [12, Lemma 14.3] in order to show the following statement formulated in Lemma B.4.
Lemma B.3
([12]) Let \(\gamma >1\), \({\underline{r}} = \frac{1}{2} \min \limits _{(t,x) \in Q_T} r >0 \) and \({\overline{r}} = 2 \max \limits _{(t,x) \in Q_T} r\). Then there exists \( C= C({\underline{r}}, {\overline{r}}) >0\) such that
where \( \mathbb {E}(\varrho |r) = P(\varrho ) - P'(r) (\varrho -r) - P(r) \) and
Now we are ready to show the following lemma.
Lemma B.4
Let \((\varrho _h, \varvec{u}_h)\) be a solution obtained either by the FV method (2.10) or the MAC method (2.11), and let \(\varvec{U}\in L^\infty ( 0,T; W^{2,\infty }( \mathbb {T}^d; \mathbb {R}^d ))\). Then there holds
where \(C_1\), \(C_2\) are the same as in Lemma B.2, and \(C_0\) depends on \({\underline{r}}, {\overline{r}}, \delta \), \( M, E_0, \gamma \).
Proof
First, thanks to Lemma B.3 we observe
where \(C=C({\underline{r}}, {\overline{r}})\) is given in Lemma B.3.
Next, using the triangular inequality, Young’s inequality, the above estimate, Lemma B.2 and Lemma B.3 we find
where \(C_0\) depends on \({\underline{r}}\), \(C({\underline{r}}, \overline{ r})\), \(\delta \), and \(C_1\). We have completed the proof. \(\square \)
Relative energy norm
In this section we show how to control the errors in the conservative variables by the relative energy.
Lemma C.1
Let \(\gamma >1\) and \((r, \varvec{U})\) satisfy
for some positive constants \(\overline{u}, {\underline{r}}, \overline{r}.\)
-
If \(\varrho >0\) and \(\int _{\mathbb {T}^d} \varrho ^\gamma \,\textrm{d}x \le E_0\) hold, then
$$\begin{aligned}{} & {} \Vert \varrho - r\Vert _{L^\gamma } + \Vert \varvec{m}- \varvec{M}\Vert _{L^{\frac{2\gamma }{\gamma +1}}} \lesssim \left( {\mathfrak {E}}(\varrho ,\varvec{u}|r,\varvec{U}) \right) ^{\frac{1}{2} } + \left( {\mathfrak {E}}(\varrho ,\varvec{u}|r,\varvec{U}) \right) ^{\frac{1}{\gamma }} \text{ for } \gamma \le 2;\nonumber \\ \end{aligned}$$(C.1a)$$\begin{aligned}{} & {} \Vert \varrho - r\Vert _{L^2} + \Vert \varvec{m}- \varvec{M}\Vert _{L^{\frac{2\gamma }{\gamma +1}}} \lesssim \left( {\mathfrak {E}}(\varrho ,\varvec{u}|r,\varvec{U}) \right) ^{\frac{1}{2} } \text{ for } \gamma \ge 2, \end{aligned}$$(C.1b)where \(m= \varrho \varvec{u}\) and \(\varvec{M}= r \varvec{U}\).
-
In addition, let \( \varrho <\overline{\varrho }\). Then
$$\begin{aligned} \Vert \varrho - r\Vert _{L^2} + \Vert \varvec{m}- \varvec{M}\Vert _{L^2} \lesssim \left( {\mathfrak {E}}(\varrho ,\varvec{u}|r,\varvec{U}) \right) ^{\frac{1}{2} }. \end{aligned}$$(C.2)
Proof
First, by the triangular inequality and Lemma B.3 we obtain for \(\gamma \le 2\) that
where \(1_{\textrm{ess}}(\varrho )\) and \(1_{\textrm{res}}(\varrho )\) are given in Lemma B.3. Further, utilizing the above estimate with the triangular inequality, Hölder’s inequality, and the \(L^\gamma \) bound on \(\varrho \), we find
which proves (C.1a).
Next, again by the triangular inequality and Lemma B.3 we observe for \(\gamma \ge 2\) that
where we have used the fact that \(\varrho ^2 \le \varrho ^\gamma \) for large \(\varrho \) with \(\gamma \ge 2\). Further, it is easy to check that
which proves (C.1b).
When assuming an upper bound on \(\varrho \), we derive via Lemma B.3 that
which implies
Combining the above two estimates we get (C.2) and complete the proof. \(\square \)
Derivation of the relative energy
In this section we show the relative energy inequality (3.5). We start with the reformulation of the relative energy.
where
Next, we collect the energy estimate (2.13), and set the test function \(\phi = \left( \frac{1}{2} \left|\varvec{U} \right|^2 - {P} '(r) \right) \) in the consistency formulation (2.15a), as well as \(\varvec{\phi }=- \varvec{U}\) in the consistency formulation (2.15b) to get respectively the following
Moreover, the term \(T_4\) reads
Summing up the above terms we get (3.5)
Rights and permissions
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.
About this article
Cite this article
Feireisl, E., Lukáčová-Medvid’ová, M. & She, B. Improved error estimates for the finite volume and the MAC schemes for the compressible Navier–Stokes system. Numer. Math. 153, 493–529 (2023). https://doi.org/10.1007/s00211-023-01346-y
Received:
Accepted:
Published:
Issue Date:
DOI: https://doi.org/10.1007/s00211-023-01346-y