Entropy Correct Spatial Discretizations for 1D Regularized Systems of Equations for Gas Mixture Dynamics

Abstract

One-dimensional regularized systems of equations for the general (multi-velocity and multi-temperature) and one-velocity and one-temperature compressible multicomponent gas mixture dynamics are considered in the absence of chemical reactions. Two types of the regularization are taken. For the latter system, diffusion fluxes between the components of the mixture are taken into account. For both the systems, the important mixture entropy balance equations with non-negative entropy productions are valid. By generalizing a discretization constructed previously in the case of a single-component gas, we suggest new nonstandard symmetric three-point spatial discretizations for both the systems which are not only conservative in mass, momentum, and total energy but also satisfy semi-discrete counterparts of the mentioned entropy balance equations with non-negative entropy productions. Importantly, the basic discretization in the one-velocity and one-temperature case is not constructed directly but by aggregation of the discretization in the case of general mixture, and that is a new approach. In this case, the results of numerical experiments are also presented for contact problems between two different gases for initial pressure jumps up to 2500.

1. Introduction

The compressible multicomponent gas mixture dynamics is of great theoretical and applied interest, and there exist various complicated systems of partial differential equations to describe it under different assumptions; in particular, see [1,2,3,4,5,6]. An important point of both physical and mathematical theories of these equations is the fulfillment of the entropy balance equation with non-negative entropy production, which confirms their physical correctness and makes it possible to prove basic a priori estimates of solutions.

A large amount of numerical methods for solving single-component gas dynamics systems of equations were developed; in particular, see [7,8,9]. Among them, there are various methods based on preliminary regularizations of these equations. These include methods based on the so-called quasi-gasdynamic (QGD) and quasi-hydrodynamic (QHD) regularizations presented in detail in monographs [10,11,12] and many subsequent papers, as well as on other regularizations [13,14,15], etc. The QGD and QHD equations for the general (multi-velocity and multi-temperature) as well as one-velocity and one-temperature gas mixture dynamics were developed, studied and applied, in particular, in [11,16,17,18,19,20,21,22], including the validity of the entropy balance equation with non-negative entropy production in [16,19,22].

The construction of entropy correct discretizations of gas dynamics systems of equations is of great interest; in particular, see [23,24,25,26,27] and references therein. In [28], a symmetric three-point version of the 1D QGD-regularized gas dynamics equations for a perfect polytropic single-component gas was constructed, which was conservative in the mass, momentum, and total energy and satisfied the semi-discrete entropy balance equation with a non-negative entropy production. This was achieved, in particular, by introducing non-standard nonlinear averages of the sought functions. In this discretization, the main sought functions (the density, velocity, and specific internal energy of the gas) were defined on a common main mesh, whereas the mass, momentum, and total energy fluxes were defined on another auxiliary mesh. Subsequently, similar discretizations were developed in the case of general equations of state of a gas [29] as well as the multidimensional case [30]. Practical testing of these discretizations was accomplished in [29,31]. Recall that it is well-known that a regularization is necessary to ensure stability of symmetric spatial discretizations combined with an explicit discretization in time, and the QGD and QHD regularizations applied by us are not trivial and physically motivated [10,11,12].

In this paper, we consider 1D regularized systems of equations for general, multi-velocity and multi-temperature, as well as one-velocity and one-temperature compressible multicomponent gas mixture dynamics of perfect polytropic gases in the absence of chemical reactions. The mixtures are supposed to be homogeneous; i.e., all the components occupy one and the same volume [2]. In the latter case, where the velocities and temperatures (but not densities) of all components are assumed to be the same, diffusion fluxes between mixture components of type [1] are taken into account. We develop the mentioned discretization from [28] further for these systems. The latter system is of primary interest for us in this paper, but first, we need to consider the former system and its discretization as well, in accordance with our aggregation approach. The systems like our main one can be applied to solve numerically gas dynamics problems with contacts between different gases [19,20,21], etc. Vast literature is devoted to problems of such kind, and we refer the reader to the recent paper [32] for a brief review and a lot of references. Note that the study of the entropy properties of numerical methods for gas mixtures is a complicated matter, and only a few recent papers deal with this subject [33].

For both the systems, we construct non-standard symmetric three-point spatial discretizations, which are not only conservative in the mass, momentum, and total energy but also satisfy the semi-discrete entropy balance equations with non-negative entropy productions. Importantly, in the one-velocity and one-temperature case, the basic discretization is not constructed directly but rather by aggregation of the discretization in the case of general mixture, and that is a new approach. Additionally, an adequate discretization is accomplished for the terms associated with the diffusion fluxes between the mixture components to guarantee the non-negative entropy production. Notice that both the QGD and QHD regularizations are considered in a unified way (by introducing a parameter into the equations like in [22]) that essentially shortens the overall study.

The paper is structured as follows. In Section 2, the 1D regularized system of equations for the dynamics of the general (multi-velocity and multi-temperature) multicomponent gas mixture is presented in the form of mass, momentum, and total energy balance equations for the mixture components. For it, the balance equations for the kinetic and internal energies of the components together with the mixture entropy balance equation with the non-negative entropy production are given. The balance equations for the mass, momentum, and total energy of the mixture are also written. In Section 3, the 1D regularized system of equations for the dynamics of the one-velocity and one-temperature multicomponent gas mixture is given, additionally taking into account diffusion fluxes between the components. A connection between the solutions to the equations for dynamics of a single-component gas and such a mixture of gases is indicated. For the one-velocity and one-temperature mixture, the corresponding balance equations for the kinetic and internal energies of the mixture as well as the mixture entropy balance equation with a non-negative entropy production are given.

In Section 4, first, meshes, mesh operators and a set of formulas for the mesh analysis are introduced. Then, a non-standard three-point spatially symmetric discretization of the regularized system of equations for the general multicomponent gas mixture dynamics is performed. Semi-discrete balance equations for the kinetic and internal energies of the components are given. Most importantly, the semi-discrete balance equation for the mixture entropy with a non-negative entropy production is derived in Theorem 1. Semi-discrete balance equations for the momentum and the total energy of the mixture are also presented. In Section 5, based on the component mass balance equation and the last two mentioned equations, a discretization of the regularized system of equations for the dynamics of the one-velocity mixture of gases is constructed including the diffusion fluxes between the components of the mixture. A connection between the solutions to the semi-discrete equations for dynamics of a single-component gas and the homogeneous mixture of gases is described. Semi-discrete balance equations for the mass, kinetic and internal energies of the mixture are given. The main result (Theorem 2) concerns the derivation of a semi-discrete balance equation for the one-velocity and one-temperature mixture entropy with a non-negative entropy production. Below, we omit words “and one-temperature” for brevity. The final Section 6 is devoted to numerical experiments in the case of the one-velocity binary mixture, more precisely, the contact between two different gases. Numerical examples concern three cases of subsonic and partially supersonic flows with shock waves, for initial data with pressure jumps up to 2500, well-known in the literature [34,35,36], and both the QGD and QHD regularizations are tested. Appendix A provides a proof of an important representation for the regularizing terms (together with heat source terms) in the semi-discrete entropy production, ensuring their non-negativity.

2. One-Dimensional (1D) Regularized System of Equations for Dynamics of General Multicomponent Gas Mixtures

A spatially 1D regularized system of equations describes the dynamics of a general multicomponent gas mixture and consists of the following balance equations for the mass, momentum, and total energy of the mixture components

$$\begin{array}{c}{\partial }_t{\rho }_{\alpha }+{\partial }_x\left({\rho }_{\alpha }(u_{\alpha }-w_{\ell \alpha })\right)=0,\end{array}$$

(1)

$$\begin{array}{c}{\partial }_t\left({\rho }_{\alpha }{u}_{\alpha }\right)+{\partial }_x\left\{{\rho }_{\alpha }(u_{\alpha }-w_{\ell \alpha })u_{\alpha }\right\}+{\partial }_x{p}_{\alpha }={\partial }_x{\mathsf{\Pi }}_{\ell \alpha }+S_{u,\alpha }+\left({\rho }_{\alpha }-\ell \tau {\partial }_x\left({\rho }_{\alpha }{u}_{\alpha }\right)\right)f_{\alpha },\end{array}$$

(2)

$$\begin{array}{c}{\partial }_t{E}_{\alpha }+{\partial }_x\left\{(E_{\alpha }+p_{\alpha })(u_{\alpha }-w_{\ell \alpha })\right\}={\partial }_x(-q_{\ell \alpha }+{\mathsf{\Pi }}_{\ell \alpha }{u}_{\alpha })+S_{E,\alpha }+{\rho }_{\alpha }(u_{\alpha }-w_{\ell \alpha })f_{\alpha }+Q_{\alpha },\end{array}$$

(3)

where $\alpha =\overline{1,K}$ (i.e., $\alpha =1,\dots ,K$) and $K⩾2$ is the amount of mixture components. The main sought functions ${\rho }_{\alpha }>0$, $u_{\alpha }$ and ${\theta }_{\alpha }>0$ are the density, velocity and absolute temperature of the mixture component $\alpha$, $\alpha =\overline{1,K}$, depending on $(x,t)$. Thus, it is assumed that each gas component has its own velocity and temperature. Here ${\partial }_t$ and ${\partial }_x$ are the partial derivatives in $t\in [0,T]$ and $x\in [-X,X]$; the cases $x\in \mathbb{R}$, $[-X,+\infty )$ or $(-\infty ,X]$ are similar. We use the perfect polytropic gas equations of state and the standard expression for the total energy of the components

$$\begin{array}{c}{p}_{\alpha }=({\gamma }_{\alpha }-1){\rho }_{\alpha }{\epsilon }_{\alpha }=R_{\alpha }{\rho }_{\alpha }{\theta }_{\alpha },{\epsilon }_{\alpha }=c_{V\alpha }{\theta }_{\alpha },E_{\alpha }=0.5{\rho }_{\alpha }{u}_{\alpha }^2+{\rho }_{\alpha }{\epsilon }_{\alpha },\end{array}$$

(4)

where $p_{\alpha }$ and ${\epsilon }_{\alpha }$ are the pressure and specific internal energy of the component $\alpha$, with the physical constants ${\gamma }_{\alpha }>1$, $c_{V\alpha }>0$ and $R_{\alpha }=({\gamma }_{\alpha }-1)c_{V\alpha }$, $\alpha =\overline{1,K}$. In addition, $R_{\alpha }=\frac{R_0}{m_{\alpha }}$, where $R_0$ is the universal gas constant (the Boltzmann constant) and $m_{\alpha }$ is the molecular mass of gas $\alpha$.

This regularization is the so-called quasi-gasdynamic (QGD) for $\ell =1$ or simpler quasi-hydrodynamic (QHD) for $\ell =0$ [16]. Hereafter, by introducing the parameter ℓ, we can treat both regularizations in a unified manner. The more universal former one is used for flows at any Mach numbers, while the latter is mainly used for subsonic and transonic flows. The regularizing velocities have the form

$$\begin{array}{c}{w}_{\ell \alpha }=\ell \frac{\tau }{{\rho }_{\alpha }}{u}_{\alpha }{\partial }_x\left({\rho }_{\alpha }{u}_{\alpha }\right)+{\widehat{w}}_{\alpha },{\widehat{w}}_{\alpha }=\frac{\tau }{{\rho }_{\alpha }}({\rho }_{\alpha }{u}_{\alpha }{\partial }_x{u}_{\alpha }+{\partial }_x{p}_{\alpha }-{\rho }_{\alpha }{f}_{\alpha }).\end{array}$$

(5)

The viscous stress ${\mathsf{\Pi }}_{\ell \alpha }$ and the heat flux $q_{\ell \alpha }$ are given by the formulas

$$\begin{array}{c}{\mathsf{\Pi }}_{\ell \alpha }={\nu }_{\alpha }{\partial }_x{u}_{\alpha }+{\mathsf{\Pi }}_{\ell \alpha }^{\tau },{\mathsf{\Pi }}_{\ell \alpha }^{\tau }={\rho }_{\alpha }{u}_{\alpha }{\widehat{w}}_{\alpha }+\ell \tau \{u_{\alpha }{\partial }_x{p}_{\alpha }+{\gamma }_{\alpha }{p}_{\alpha }{\partial }_x{u}_{\alpha }-({\gamma }_{\alpha }-1)Q_{\alpha }\},\end{array}$$

(6)

$$\begin{array}{c}{q}_{\ell \alpha }=-{\varkappa }_{\alpha }{\partial }_x{\theta }_{\alpha }+\ell q_{\alpha }^{\tau },-q_{\alpha }^{\tau }=\tau \left\{{\rho }_{\alpha }{u}_{\alpha }^2\left({\partial }_x{\epsilon }_{\alpha }-\frac{p_{\alpha }}{{\rho }_{\alpha }^2}{\partial }_x{\rho }_{\alpha }\right)-Q_{\alpha }{u}_{\alpha }\right\},\end{array}$$

(7)

where ${\nu }_{\alpha }⩾0$, ${\varkappa }_{\alpha }⩾0$ and $\tau >0$ are the coefficients of viscosity and heat conductivity (physical or artificial) and the relaxation (regularization) parameter which may depend on the sought functions (here, their specific form is not essential). The terms ${\nu }_{\alpha }{\partial }_x{u}_{\alpha }$ and $-{\varkappa }_{\alpha }{\partial }_x{\theta }_{\alpha }$ are of the Navier–Stokes–Fourier type, whereas the terms ${\mathsf{\Pi }}_{\ell \alpha }^{\tau }$ and $q_{\alpha }^{\tau }$ are relaxation (regularizing) ones. The densities of body forces $f_{\alpha }$ and intensities of heat sources $Q_{\alpha }⩾0$ are given functions.

In the case of the binary mixture ($K=2$), such equations were considered in [16] (in the multidimensional version). The above equations for the multicomponent mixture generalize them in a natural way. Note that

$$\begin{array}{c}{w}_{1\alpha }=\frac{\tau }{{\rho }_{\alpha }}\left({\partial }_x\left({\rho }_{\alpha }{u}_{\alpha }^2\right)+{\partial }_x{p}_{\alpha }-{\rho }_{\alpha }{f}_{\alpha }\right),w_{0\alpha }={\widehat{w}}_{\alpha }.\end{array}$$

(8)

The written equations form a set of standard regularized equations for the balance of mass, momentum, and total energy of individual components of the mixture, where $S_{u,\alpha }$ and $S_{E,\alpha }$ are the additional exchange terms. They depend on the sought functions ${\{{\rho }_{\alpha },u_{\alpha },{\theta }_{\alpha }\}}_{\alpha =1}^K$ and satisfy the relations

$$\begin{array}{c}\langle S_{u,\alpha }\rangle =0,\langle S_{u,\alpha }{u}_{\alpha }\rangle ⩽0,\langle S_{E,\alpha }\rangle =0,〈\frac{1}{{\theta }_{\alpha }}(S_{E,\alpha }-S_{u,\alpha }{u}_{\alpha })〉⩾0.\end{array}$$

(9)

Hereafter, $\langle ·\rangle$ denotes the summation operation over the index $\alpha =\overline{1,K}$. The second relation is not used in this paper.

In the case of the binary mixtue ($K=2$), for any ${\gamma }_{\alpha }>1$, $\alpha =1,2$, the exchange terms can be specified, in particular, by means of formulas of the molecular–kinetic type [16]

$$\begin{array}{c}{S}_{u,\alpha }={\nu }_{\alpha \beta }{\rho }_{\alpha }(\widehat{u}-u_{\alpha }),S_{E,\alpha }={\nu }_{\alpha \beta }({\widehat{E}}_{\alpha }-E_{\alpha }),\widehat{u}=\frac{m_1{u}_1+m_2{u}_2}{m_1+m_2},{\widehat{E}}_{\alpha }=\frac{1}{2}{\rho }_{\alpha }{\widehat{u}}^2+c_{V\alpha }{\rho }_{\alpha }{\widehat{\theta }}_{\alpha },\\ {\widehat{\theta }}_{\alpha }={\theta }_{\alpha }+C^{\left(m\right)}({\gamma }_{\alpha }-1)\left\{\frac{3}{2}({\theta }_{\beta }-{\theta }_{\alpha })+\frac{m_{\beta }}{4R_0}{(u_1-u_2)}^2\right\},C^{\left(m\right)}=\frac{2m_1{m}_2}{{(m_1+m_2)}^2},\end{array}$$

where ${\nu }_{\alpha \beta }>0$ is the frequency with which molecules of gas $\alpha$ collide (bump into) molecules of gas $\beta =3-\alpha$, $\alpha =1,2$. Here, it is assumed that the following balance equality holds

$$\frac{{\rho }_1}{m_1}{\nu }_{12}=\frac{{\rho }_2}{m_2}{\nu }_{21}=N_{col},$$

where $N_{col}$ is the total number of collisions between the molecules of gases 1 and 2, and consequently, relations (9) are valid; see their derivation in [16]. Formulas for ${\nu }_{12}$ and ${\nu }_{21}$ can be found in [11,37] and are not reproduced here.

In the general case $K⩾2$, in particular, the formulas of type [38] can be used

$$S_{u,\alpha }={\langle K_{\alpha \beta }(u_{\beta }-u_{\alpha })\rangle }_{\beta },S_{E,\alpha }={\langle c_{\alpha \beta }({\theta }_{\beta }-{\theta }_{\alpha })\rangle }_{\beta }+{\langle b_{\alpha \beta }{K}_{\alpha \beta }{(u_{\beta }-u_{\alpha })}^2\rangle }_{\beta }+S_{u,\alpha }{u}_{\alpha }.$$

Here, ${\langle ·\rangle }_{\beta }$ is the summation operation over the index $\beta =\overline{1,K}$.

For $K_{\alpha \beta }=K_{\beta \alpha }$, $c_{\alpha \beta }=c_{\beta \alpha }$ and $b_{\alpha \beta }+b_{\beta \alpha }=1$ for any $\alpha \ne \beta$, swapping the indexes $\alpha$ and $\beta$, one derives

$$\begin{array}{c}\langle {\langle K_{\alpha \beta }(u_{\beta }-u_{\alpha })\rangle }_{\beta }\rangle =\langle {\langle K_{\alpha \beta }(u_{\alpha }-u_{\beta })\rangle }_{\beta }\rangle =0,\langle {\langle c_{\alpha \beta }({\theta }_{\beta }-{\theta }_{\alpha })\rangle }_{\beta }\rangle =0,\\ \langle {\langle K_{\alpha \beta }(u_{\beta }-u_{\alpha })\rangle }_{\beta }{u}_{\alpha }\rangle =\frac{1}{2}\langle \langle {K_{\alpha \beta }\left((u_{\beta }-u_{\alpha })u_{\alpha }+(u_{\alpha }-u_{\beta })u_{\beta }\right){\rangle }_{\beta }\rangle =-\frac{1}{2}\langle \langle K_{\alpha \beta }{(u_{\beta }-u_{\alpha })}^2\rangle }_{\beta }\rangle ,\\ \langle \langle {b_{\alpha \beta }{K}_{\alpha \beta }{(u_{\beta }-u_{\alpha })}^2\rangle }_{\beta }\rangle =\frac{1}{2}\langle \langle {(b_{\alpha \beta }+b_{\beta \alpha })K_{\alpha \beta }{(u_{\beta }-u_{\alpha })}^2\rangle }_{\beta }\rangle =\frac{1}{2}\langle \langle {K_{\alpha \beta }{(u_{\beta }-u_{\alpha })}^2\rangle }_{\beta }\rangle ,\\ \langle \frac{1}{{\theta }_{\alpha }}{\langle c_{\alpha \beta }({\theta }_{\beta }-{\theta }_{\alpha })\rangle }_{\beta }\rangle =\frac{1}{2}\langle {\langle c_{\alpha \beta }(\frac{1}{{\theta }_{\alpha }}({\theta }_{\beta }-{\theta }_{\alpha })+\frac{1}{{\theta }_{\beta }}({\theta }_{\alpha }-{\theta }_{\beta }))\rangle }_{\beta }\rangle \\ =\frac{1}{2}\langle {\langle \frac{c_{\alpha \beta }}{{\theta }_{\alpha }{\theta }_{\beta }}{({\theta }_{\beta }-{\theta }_{\alpha })}^2\rangle }_{\beta }\rangle .\end{array}$$

Consequently, relations (9) are valid once again provided that $K_{\alpha \beta }⩾0$, $b_{\alpha \beta }⩾0$ and $c_{\alpha \beta }⩾0$ for any $\alpha \ne \beta$.

In the absence of the regularization, i.e., for $\tau =0$, the above regularized system of equations is reduced to the 1D compressible Navier–Stokes–Fourier-type system for multicomponent flows for ${\nu }_{\alpha }>0$ and ${\varkappa }_{\alpha }>0$ or the Euler-type one for ${\nu }_{\alpha }={\varkappa }_{\alpha }=0$, $\alpha =\overline{1,K}$; in particular, see [6,38,39,40] and references therein.

Now, we present some important balance equations that follow from the basic ones.

Lemma 1.

Let $j_{\ell \alpha }={\rho }_{\alpha }(u_{\alpha }-w_{\ell \alpha })$ be the regularized mass flux of the mixture component α. The following balance equations for the kinetic and internal energies of the mixture components hold

$$\begin{array}{c}0.5{\partial }_t\left({\rho }_{\alpha }{u}_{\alpha }^2\right)+0.5{\partial }_x\left(j_{\ell \alpha }{u}_{\alpha }^2\right)+\left({\partial }_x{p}_{\alpha }\right)u_{\alpha }\\ =\left({\partial }_x{\mathsf{\Pi }}_{\ell \alpha }\right)u_{\alpha }+S_{u,\alpha }{u}_{\alpha }+\left({\rho }_{\alpha }-\ell \tau {\partial }_x\left({\rho }_{\alpha }{u}_{\alpha }\right)\right)f_{\alpha }{u}_{\alpha },\end{array}$$

(10)

$$\begin{array}{c}{\partial }_t\left({\rho }_{\alpha }{\epsilon }_{\alpha }\right)+{\partial }_x\left(j_{\ell \alpha }{\epsilon }_{\alpha }\right)+p_{\alpha }{\partial }_x{u}_{\alpha }={\partial }_x(-q_{\ell \alpha }+p_{\alpha }{w}_{\ell \alpha })\\ +{\mathsf{\Pi }}_{\ell \alpha }{\partial }_x{u}_{\alpha }+S_{E,\alpha }-S_{u,\alpha }{u}_{\alpha }-{\rho }_{\alpha }{\widehat{w}}_{\alpha }{f}_{\alpha }+Q_{\alpha }\end{array}$$

(11)

for $(x,t)\in [-X,X]\times [0,T]$, where $\alpha =\overline{1,K}$.

Proof. 

To derive these standard type equations, we first multiply the momentum balance equation for the components (2) by $u_{\alpha }$ and use the formulas of the same type

$$\begin{array}{c}\left({\partial }_t\left({\rho }_{\alpha }{u}_{\alpha }\right)\right)u_{\alpha }=0.5{\partial }_t\left({\rho }_{\alpha }{u}_{\alpha }^2\right)+0.5\left({\partial }_t{\rho }_{\alpha }\right)u_{\alpha }^2,\\ \left({\partial }_x\left(j_{\alpha }{u}_{\alpha }\right)\right)u_{\alpha }=0.5{\partial }_x\left(j_{\alpha }{u}_{\alpha }^2\right)+0.5\left({\partial }_x{j}_{\alpha }\right)u_{\alpha }^2.\end{array}$$

(12)

Applying the mass balance equation for the components (1), we obtain Equation (10).

We subtract it from the balance equation for the total energy of the components (3). By differentiating with respect to x the products $p_{\alpha }{u}_{\alpha }$ and ${\mathsf{\Pi }}_{\ell \alpha }{u}_{\alpha }$ and using the first Formula (5), we derive the balance equation for the internal energy of the components (11). □

Let us introduce the total density, the specific entropies of the components and the specific total entropy (the mixture entropy)

$$\rho =\langle {\rho }_{\alpha }\rangle ,s_{\alpha }={\overline{s}}_{\alpha }-R_{\alpha }ln\frac{{\rho }_{\alpha }}{{\overline{\rho }}_{\alpha }}+c_{V\alpha }ln\frac{{\theta }_{\alpha }}{{\overline{\theta }}_{\alpha }},s=〈\frac{{\rho }_{\alpha }}{\rho }{s}_{\alpha }〉,$$

(13)

where ${\overline{s}}_{\alpha }$, ${\overline{\rho }}_{\alpha }>0$ and ${\overline{\theta }}_{\alpha }>0$ are constants (reference values of $s_{\alpha }$, ${\rho }_{\alpha }$ and ${\theta }_{\alpha }$). Here, $C_{\alpha }:=\frac{{\rho }_{\alpha }}{\rho }$ are the mass concentrations of the mixture components; they will arise below as well.

The entropy balance equation for the general gas mixture has the form

$$\begin{array}{c}{\partial }_t\left(\rho s\right)+{\partial }_x\langle j_{\ell \alpha }{s}_{\alpha }\rangle ={\partial }_x〈\frac{1}{{\theta }_{\alpha }}({\varkappa }_{\alpha }{\partial }_x{\theta }_{\alpha }-\ell q_{\alpha }^{\tau })〉\\ +〈\frac{{\varkappa }_{\alpha }}{{\theta }_{\alpha }^2}{\left({\partial }_x{\theta }_{\alpha }\right)}^2+\frac{{\nu }_{\alpha }}{{\theta }_{\alpha }}{\left({\partial }_x{u}_{\alpha }\right)}^2〉+〈\frac{1}{{\theta }_{\alpha }}(S_{E,\alpha }-S_{u,\alpha }{u}_{\alpha })〉\\ +〈\frac{{\rho }_{\alpha }}{\tau {\theta }_{\alpha }}{\widehat{w}}_{\alpha }^2+\ell \tau \frac{R_{\alpha }}{{\rho }_{\alpha }}{\left\{{\partial }_x\left({\rho }_{\alpha }{u}_{\alpha }\right)\right\}}^2+\ell \tau \frac{c_{V_{\alpha }}{\rho }_{\alpha }}{{\theta }_{\alpha }^2}{\left\{({\gamma }_{\alpha }-1){\theta }_{\alpha }{\partial }_x{u}_{\alpha }+u_{\alpha }{\partial }_x{\theta }_{\alpha }-\frac{Q_{\alpha }}{2c_{V\alpha }{\rho }_{\alpha }}\right\}}^2〉\\ +〈\frac{Q_{\alpha }}{{\theta }_{\alpha }}\left(1-\ell \tau \frac{({\gamma }_{\alpha }-1)Q_{\alpha }}{4p_{\alpha }}\right)〉\end{array}$$

(14)

for $(x,t)\in [-X,X]\times [0,T]$. For $K=2$, it was derived in [16] (in the multidimensional case), and the proof was based on a similar equation in the single-component case ($K=1$) [11,12] (see also [41]); the derivation in the case $K⩾2$ is actually the same.

The sum of all terms on the right-hand side of Equation (14), except for the divergent first one, is the entropy production (or the dissipative function). Among these terms, the sum

$$〈\frac{{\varkappa }_{\alpha }}{{\theta }_{\alpha }^2}{\left({\partial }_x{\theta }_{\alpha }\right)}^2+\frac{{\nu }_{\alpha }}{{\theta }_{\alpha }}{\left({\partial }_x{u}_{\alpha }\right)}^2〉+〈\frac{1}{{\theta }_{\alpha }}(S_{E,\alpha }-S_{u,\alpha }{u}_{\alpha })〉+〈\frac{Q_{\alpha }}{{\theta }_{\alpha }}〉$$

is the Navier–Stokes–Fourier entropy production, and all its three summands (each in the angular brackets) are non-negative, taking into account the last property (9).

The remaining terms with the factor $\tau$ are relaxation ones. The penultimate term on the right in Equation (14) contains the sum of three quadratic terms with positive coefficients under the sign $\langle ·\rangle$, so that term is non-negative. The last term on the right in Equation (14) is non-negative for $\ell =0$, as well as for $\ell =1$ and under the conditions $\tau ({\gamma }_{\alpha }-1)Q_{\alpha }⩽4p_{\alpha }$, $\alpha =\overline{1,K}$, and then, the total entropy production is also non-negative.

This entropy balance equation and the indicated property of the entropy production non-negativity remain valid for $\tau ⩾0$, where one should switch to another form for the first relaxation term:

$$\frac{{\rho }_{\alpha }}{\tau {\theta }_{\alpha }}{\widehat{w}}_{\alpha }^2=\frac{\tau }{{\rho }_{\alpha }{\theta }_{\alpha }}{({\rho }_{\alpha }u{\partial }_{x}u+{\partial }_x{p}_{\alpha }-{\rho }_{\alpha }f)}^2.$$

Recall that with the help of the specific entropy of the components, the expression for the heat flux (7) can be written in a shortened form

$$-q_{\ell \alpha }={\varkappa }_{\alpha }{\partial }_x{\theta }_{\alpha }+\tau ({\rho }_{\alpha }{\theta }_{\alpha }{u}_{\alpha }^2{\partial }_x{s}_{\alpha }-u_{\alpha }{Q}_{\alpha });$$

however, this form is not applied below.

Applying the operation $\langle ·\rangle$ to the balance Equations (1)–(3) for the mass, momentum, and total energy of the mixture components (i.e., summing them over the index $\alpha =\overline{1,K}$), we obtain their corollaries which are the balance equations for the mass, momentum and total energy of the mixture

$$\begin{array}{c}{\partial }_t\rho +{\partial }_x\langle {\rho }_{\alpha }(u_{\alpha }-w_{\ell \alpha })\rangle =0,\end{array}$$

(15)

$$\begin{array}{c}{\partial }_t\langle {\rho }_{\alpha }{u}_{\alpha }\rangle +{\partial }_x\langle {\rho }_{\alpha }(u_{\alpha }-w_{\ell \alpha })u_{\alpha }\rangle +{\partial }_{x}p={\partial }_x{\mathsf{\Pi }}_{\ell }+〈\left({\rho }_{\alpha }-\ell \tau {\partial }_x\left({\rho }_{\alpha }{u}_{\alpha }\right)\right)f_{\alpha }〉,\end{array}$$

(16)

$$\begin{array}{c}{\partial }_{t}E+{\partial }_x\langle (E_{\alpha }+p_{\alpha })(u_{\alpha }-w_{\ell \alpha })\rangle ={\partial }_x(-q_{\ell }+\langle {\mathsf{\Pi }}_{\ell \alpha }{u}_{\alpha }\rangle )+\langle {\rho }_{\alpha }(u_{\alpha }-w_{\ell \alpha })f_{\alpha }\rangle +\langle Q_{\alpha }\rangle .\end{array}$$

(17)

Herein, the aggregated (or total, i.e., for the mixture) pressure, total energy, viscous stress and heat flux are expressed by the formulas

$$\begin{array}{c}p=\langle p_{\alpha }\rangle ,E=\langle E_{\alpha }\rangle ,{\mathsf{\Pi }}_{\ell }=\langle {\mathsf{\Pi }}_{\ell \alpha }\rangle ,q_{\ell }=\langle q_{\ell \alpha }\rangle .\end{array}$$

(18)

The first of these is Dalton’s law for mixtures. The exchange terms have been cancelled out due to the equalities $\langle S_{u,\alpha }\rangle =0$ and $\langle S_{E,\alpha }\rangle =0$, see (9). These equations are used in the next section.

Recall that, for the purposes of discretization, Equations (1) are often replaced by Equation (15) and equations for $K-1$ concentrations

$${\partial }_t\left(\rho C_{\alpha }\right)+{\partial }_x\left(\rho C_{\alpha }(u_{\alpha }-w_{\ell \alpha })\right)=0,\alpha =\overline{1,K-1},$$

although in this paper, we will not use this approach.

3. One-Dimensional (1D) Regularized System of Equations for the Dynamics of One-Velocity Multicomponent Gas Mixtures in the Presence of Diffusion Fluxes

In this section, we consider the model of the one-velocity mixture, where the velocities and temperatures of all components are assumed to be the same, which drastically reduces the number of sought functions. Let also $f_{\alpha }=f$, $\alpha =\overline{1,K}$. Additionally, diffusion fluxes between the mixture components are introduced.

As in [19,22], we write the aggregated regularized system of equations for the one-velocity mixture dynamics (in the diffusion approximation), which is obtained from Equations (1), (16) and (17) after taking $u_{\alpha }=u$ and ${\theta }_{\alpha }=\theta$ for $\alpha =\overline{1,K}$, and it consists of the balance equations for the mass of the component and the momentum and total energy of the mixture

$$\begin{array}{c}{\partial }_t{\rho }_{\alpha }+{\partial }_x\left({\rho }_{\alpha }(u-w_{\ell \alpha })+d_{\alpha }\right)=0,\alpha =\overline{1,K},\end{array}$$

(19)

$$\begin{array}{c}{\partial }_t\left(\rho u\right)+{\partial }_x\langle {\rho }_{\alpha }(u-w_{\ell \alpha })u\rangle +{\partial }_{x}p={\partial }_x{\mathsf{\Pi }}_{\ell }+\left(\rho -\ell \tau {\partial }_x\left(\rho u\right)\right)f,\end{array}$$

(20)

$$\begin{array}{c}{\partial }_{t}E+{\partial }_x\langle (E_{\alpha }+p_{\alpha })(u-w_{\ell \alpha })\rangle ={\partial }_x(-q_{\ell }+{\mathsf{\Pi }}_{\ell }u)+\langle {\rho }_{\alpha }(u-w_{\ell \alpha })\rangle f+Q,\end{array}$$

(21)

with the additional inclusion of the diffusion flux $d_{\alpha }$ between the component $\alpha$ and the rest of the components. Now, the main sought functions are the component densities ${\rho }_{\alpha }>0$, $\alpha =\overline{1,K}$, and their common velocity u and temperature $\theta >0$. Still, $(x,t)\in [-X,X]\times [0,T]$; also $Q=\langle Q_{\alpha }\rangle$.

The pressure, specific internal and total energies of the components take a slightly simplified form

$$p_{\alpha }=({\gamma }_{\alpha }-1){\rho }_{\alpha }{\epsilon }_{\alpha }=R_{\alpha }{\rho }_{\alpha }\theta ,{\epsilon }_{\alpha }=c_{V\alpha }\theta ,E_{\alpha }=0.5{\rho }_{\alpha }{u}^2+{\rho }_{\alpha }{\epsilon }_{\alpha },$$

(22)

and the total density and corresponding quantities appearing in the equations for the one-velocity mixture are as follows

$$\rho =\langle {\rho }_{\alpha }\rangle ,p=\langle p_{\alpha }\rangle =R\rho \theta ,\epsilon =〈\frac{{\rho }_{\alpha }}{\rho }{\epsilon }_{\alpha }〉=c_V\theta ,E=\langle E_{\alpha }\rangle =0.5\rho u^2+\rho \epsilon ,$$

(23)

with the coefficients

$$R:=〈\frac{{\rho }_{\alpha }}{\rho }{R}_{\alpha }〉,c_V:=〈\frac{{\rho }_{\alpha }}{\rho }{c}_{V\alpha }〉.$$

(24)

We emphasize that here, R and $c_V$ are functions rather than constants in contrast to the single-component case (provided that there are different coefficients among $R_{\alpha }$ and $c_{V\alpha }$). The expression for p can be rewritten in the standard form

$$p=(\gamma -1)\rho \epsilon \mathrm{with}\gamma :=\frac{R}{c_V}+1.$$

Note that if ${\gamma }_1=\dots ={\gamma }_K$, then also $\gamma ={\gamma }_1$ even for non-coinciding $R_{\alpha }$.

In these equations, according to the performed aggregation procedure, we have

$$\begin{array}{c}{w}_{\ell \alpha }=\ell \frac{\tau }{{\rho }_{\alpha }}u{\partial }_x\left({\rho }_{\alpha }u\right)+{\widehat{w}}_{\alpha },{\widehat{w}}_{\alpha }=\frac{\tau }{{\rho }_{\alpha }}({\rho }_{\alpha }u{\partial }_{x}u+{\partial }_x{p}_{\alpha }-{\rho }_{\alpha }f),\end{array}$$

(25)

$$\begin{array}{c}{\mathsf{\Pi }}_{\ell }=\nu {\partial }_{x}u+{\mathsf{\Pi }}_{\ell }^{\tau },{\mathsf{\Pi }}_{\ell }^{\tau }=u\langle {\rho }_{\alpha }{\widehat{w}}_{\alpha }\rangle +\ell \tau \{u{\partial }_{x}p+\langle {\gamma }_{\alpha }{p}_{\alpha }\rangle {\partial }_{x}u-\langle ({\gamma }_{\alpha }-1)Q_{\alpha }\rangle \},\end{array}$$

(26)

$$\begin{array}{c}{q}_{\ell }=-\varkappa {\partial }_x\theta +\ell q^{\tau }+q^d,\end{array}$$

(27)

$$\begin{array}{c}-q^{\tau }=\tau \left\{u^2\left(\langle c_{V\alpha }{\rho }_{\alpha }\rangle {\partial }_x\theta -\theta {\partial }_x\langle R_{\alpha }{\rho }_{\alpha }\rangle \right)-Qu\right\}=\tau \left\{u^2\left(c_V\rho {\partial }_x\theta -\theta {\partial }_x\left(R\rho \right)\right)-Qu\right\},\end{array}$$

(28)

where $\nu :=\langle {\nu }_{\alpha }\rangle$ and $\varkappa :=\langle {\varkappa }_{\alpha }\rangle$, and the appearance of the additional term $q^d$ is related to the introduction of diffusion fluxes. We emphasize that a number of new functions (including $w_{\ell \alpha }$ and ${\widehat{w}}_{\alpha }$) here retain their notation from the previous section.

As in [19,22], the introduction of total regularizing velocities

$$\begin{array}{c}{w}_{\ell }:=〈\frac{{\rho }_{\alpha }}{\rho }{w}_{\ell \alpha }〉=\ell \frac{\tau }{\rho }u{\partial }_x\left(\rho u\right)+\widehat{w},\widehat{w}:=〈\frac{{\rho }_{\alpha }}{\rho }{\widehat{w}}_{\alpha }〉=\frac{\tau }{\rho }(\rho u{\partial }_{x}u+{\partial }_{x}p-\rho f)\end{array}$$

(29)

makes it possible to simplify the form of the balance equations for the momentum and total energy of the one-velocity mixture (20) and (21):

$$\begin{array}{c}{\partial }_t\left(\rho u\right)+{\partial }_x\left(\rho (u-w_{\ell })u\right)+{\partial }_{x}p={\partial }_x{\mathsf{\Pi }}_{\ell }+\left(\rho -\ell \tau {\partial }_x\left(\rho u\right)\right)f,\end{array}$$

(30)

$$\begin{array}{c}{\partial }_{t}E+{\partial }_x\left(0.5\rho u^2(u-w_{\ell })+\langle {\rho }_{\alpha }{c}_{p\alpha }\theta (u-w_{\ell \alpha })\rangle \right)\\ ={\partial }_x(-q_{\ell }+{\mathsf{\Pi }}_{\ell }u)+\rho (u-w_{\ell })f+Q,\end{array}$$

(31)

where $c_{p\alpha }=R_{\alpha }+c_{V\alpha }={\gamma }_{\alpha }{c}_{V\alpha }$. However, for the purpose of discretization, we prefer to use the original form of these equations below. Recall that $c_{V\alpha }$ and $c_{p\alpha }$ are the specific heat capacities at constant volume and pressure and $c_{p\alpha }\theta$ is the specific enthalpy of the component $\alpha =\overline{1,K}$. Moreover, in (26), we have $\langle {\rho }_{\alpha }{\widehat{w}}_{\alpha }\rangle =\rho \widehat{w}$.

In the absence of the regularization, i.e., for $\tau =0$, the presented regularized system of equations for the one-velocity mixture dynamics is reduced to the 1D compressible Navier–Stokes–Fourier-type system for one-velocity and one-temperature multicomponent flows for ${\nu }_{\alpha }>0$ and ${\varkappa }_{\alpha }>0$ or the Euler-type one for ${\nu }_{\alpha }={\varkappa }_{\alpha }=0$, $\alpha =\overline{1,K}$, in particular, see [1,4,33] and references therein.

In these equations, we define the terms related to the diffusion fluxes between the components by the following formulas

$$\begin{array}{c}-d_{\alpha }:=d_0\left({\langle {\partial }_x(G_{\alpha }-G_{\beta })\rangle }_{\beta }+b_{\alpha }{\partial }_x\theta \right)=d_0\left({\partial }_x(K{G}_{\alpha }-G)+b_{\alpha }{\partial }_x\theta \right)\mathrm{with}G:=\langle G_{\alpha }\rangle ,\end{array}$$

(32)

$$\begin{array}{c}{q}^d=\langle (G_{\alpha }+K^{-1}{b}_{\alpha }\theta )d_{\alpha }\rangle ,\end{array}$$

(33)

$$\begin{array}{c}{G}_{\alpha }:={\epsilon }_{\alpha }-s_{\alpha }\theta +\frac{p_{\alpha }}{{\rho }_{\alpha }}=(c_{p\alpha }-s_{\alpha })\theta ,s_{\alpha }={\overline{s}}_{\alpha }-R_{\alpha }ln\frac{{\rho }_{\alpha }}{{\overline{\rho }}_{\alpha }}+c_{V\alpha }ln\frac{\theta }{\overline{\theta }},\end{array}$$

(34)

where $G_{\alpha }$ and $s_{\alpha }$ are the Gibbs potential and specific entropy of the component $\alpha =\overline{1,K}$ (for example, see [42]).

The quantities $d_0⩾0$ and $b_{\alpha }$ are not specified in this paper; they can depend on the sought functions, and it is only assumed that $\langle b_{\alpha }\rangle =0$, whereas ${\overline{s}}_{\alpha }$, ${\overline{\rho }}_{\alpha }>0$ and $\overline{\theta }>0$ are constants (reference values of $s_{\alpha }$, ${\rho }_{\alpha }$ and $\theta$). The property $\langle d_{\alpha }\rangle =0$ plays an important role and immediately follows from definitions (32) and the condition $\langle b_{\alpha }\rangle$=0.

Note that the formula $-q^{\tau }=\tau \left(\theta u^2\langle {\rho }_{\alpha }{\partial }_x{s}_{\alpha }\rangle -Qu\right)$ is valid.

In the case of the binary mixture ($K=2$) with $d_{\alpha }=0$, $\alpha =1,2$, the above equations for the dynamics of the one-velocity mixture were obtained in [19] (in the multidimensional case). The general multicomponent case ($K⩾2$) with the introduction of $d_{\alpha }$ and $q^d$ has recently been considered in [22] for $\ell =0,1$. The above quantities $d_{\alpha }$ and $q^d$ generalize those proposed in [1] in the case of the binary mixture; this particular case is discussed in more detail in [22].

Now, we pass to several properties of the above system of equations for the one-velocity gas mixture dynamics.

Proposition 1.

Let the functions $\rho >0$, u and $\theta >0$ be a solution to the following regularized system of equations for a single-component gas dynamics

$$\begin{array}{c}{\partial }_t\rho +{\partial }_x\left(\rho (u-w_{\ell })\right)=0,\end{array}$$

(35)

$$\begin{array}{c}{\partial }_t\left(\rho u\right)+{\partial }_x\left(\rho (u-w_{\ell })u\right)+{\partial }_{x}p={\partial }_x{\mathsf{\Pi }}_{\ell }+\left(\rho -\ell \tau {\partial }_x\left(\rho u\right)\right)f,\end{array}$$

(36)

$$\begin{array}{c}{\partial }_{t}E+{\partial }_x\left((E+p)(u-w_{\ell })\right)={\partial }_x(-q_{\ell }+{\mathsf{\Pi }}_{\ell }u)+\rho (u-w_{\ell })f+Q\end{array}$$

(37)

for $(x,t)\in [-X,X]\times [0,T]$, where $p=(\gamma -1)\rho \epsilon =R\rho \theta$, $\epsilon =c_V\theta$ and $E=0.5\rho u^2+\rho \epsilon$ with constant $R>0$ and $c_V>0$ together with

$$\begin{array}{c}{w}_{\ell }=\ell \frac{\tau }{\rho }u{\partial }_x\left(\rho u\right)+\widehat{w},\widehat{w}=\frac{\tau }{\rho }(\rho u{\partial }_{x}u+{\partial }_{x}p-\rho f),\end{array}$$

(38)

$$\begin{array}{c}{\mathsf{\Pi }}_{\ell }=\nu {\partial }_{x}u+{\mathsf{\Pi }}_{\ell }^{\tau },{\mathsf{\Pi }}_{\ell }^{\tau }=u\rho \widehat{w}+\ell \tau \{u{\partial }_{x}p+\gamma p{\partial }_{x}u-(\gamma -1)Q\},\end{array}$$

(39)

$$\begin{array}{c}{q}_{\ell }=-\varkappa {\partial }_x\theta +\ell q^{\tau },-q^{\tau }=\tau \left\{u^2\left(c_V\rho {\partial }_x\theta -R\theta {\partial }_x\rho \right)-Qu\right\}\end{array}$$

(40)

with $\tau =\tau (\rho ,\epsilon ,u)>0$, $\nu =\nu (\rho ,\epsilon ,u)⩾0$ and $\varkappa =\varkappa (\rho ,\epsilon ,u)⩾0$.

Then, for any constant $0<C_{\alpha }<1$, $\alpha =\overline{1,K}$, such that $\langle C_{\alpha }\rangle =1$, the functions ${\rho }_{\alpha }=C_{\alpha }\rho$, $\alpha =\overline{1,K}$, u and θ are a solution to the regularized system of equations for the one-velocity mixture dynamics (19)–(28) with $d_{\alpha }=0$, ${\gamma }_{\alpha }=\gamma$, $R_{\alpha }=R$ and $c_{V\alpha }=c_V$, $\alpha =\overline{1,K}$, as well as $\tau =\tau (\rho ,\epsilon ,u)$, $\nu =\nu (\rho ,\epsilon ,u)$ and $\varkappa =\varkappa (\rho ,\epsilon ,u)$ with ρ and ε given by formulas (23).

Proof. 

Under the hypotheses of the proposition, Formulas (22)–(24) are valid. Moreover, we have ${\widehat{w}}_{\alpha }=\widehat{w}$ and $w_{\ell \alpha }=w_{\ell }$ in definition (25), see definition (38), so Formulas (39) and (40) coincide with the corresponding formulas for the mixture (26)–(28). This means that the momentum and total energy balance equations (36) and (37) coincide with the corresponding equations for the mixture (30) and (31). Finally, multiplying the mass balance Equation (35) by $C_{\alpha }$ yields the corresponding equation for the mixture components (19). □

This proposition expresses the natural fact that the formal decomposition of a single-component gas into K components with proportional densities and the same temperatures, velocities and physical constants is a special (trivial) case of the one-velocity mixture. In particular, Proposition 1 is useful for checking various properties of solutions to the equations of the one-velocity mixture dynamics. Of course, its multidimensional version holds as well.

Applying the operation $\langle ·\rangle$ to the mass balance equation for the mixture components (19) and using the property $\langle d_{\alpha }\rangle =0$, we obtain the important balance equation for the total mass

$${\partial }_t\rho +{\partial }_x\langle j_{\ell \alpha }\rangle =0,$$

(41)

where $j_{\ell \alpha }={\rho }_{\alpha }(u-w_{\ell \alpha })$ is the regularized mass flux of the component $\alpha$; moreover, $\langle j_{\ell \alpha }\rangle =\rho (u-w_{\ell })$ according to (29).

The above balance Equations (30), (41) and (31) entail the balance equations for the kinetic and internal energies of the one-velocity mixture [22]

$$\begin{array}{c}0.5{\partial }_t\left(\rho u^2\right)+0.5{\partial }_x\left(\rho (u-w_{\ell })u^2\right)+\left({\partial }_{x}p\right)u=\left({\partial }_x{\mathsf{\Pi }}_{\ell }\right)u+\left(\rho -\ell \tau {\partial }_x\left(\rho u\right)\right)fu,\end{array}$$

(42)

$$\begin{array}{c}{\partial }_t\left(\rho \epsilon \right)+{\partial }_x\langle j_{\ell \alpha }{\epsilon }_{\alpha }\rangle +p{\partial }_{x}u={\partial }_x(-q_{\ell }+\langle p_{\alpha }{w}_{\ell \alpha }\rangle )+{\mathsf{\Pi }}_{\ell }{\partial }_{x}u-\rho \widehat{w}f+Q.\end{array}$$

(43)

The specific entropy of the one-velocity mixture is still given by the formula $s=〈\frac{{\rho }_{\alpha }}{\rho }{s}_{\alpha }〉$. The entropy balance equation for the one-velocity mixture has the form [22]

$$\begin{array}{c}{\partial }_t\left(\rho s\right)+{\partial }_x\langle j_{\ell \alpha }{s}_{\alpha }\rangle ={\partial }_x\left(\frac{1}{\theta }(\varkappa {\partial }_x\theta -\ell q^{\tau })-\frac{1}{K}\langle b_{\alpha }{d}_{\alpha }\rangle \right)\\ +\frac{\varkappa }{{\theta }^2}{\left({\partial }_x\theta \right)}^2+\frac{\nu }{\theta }{\left({\partial }_{x}u\right)}^2+\frac{1}{K{d}_0\theta }\langle d_{\alpha }^2\rangle +\frac{1}{\tau \theta }\langle {\rho }_{\alpha }{\widehat{w}}_{\alpha }^2\rangle +\ell \tau 〈\frac{R_{\alpha }}{{\rho }_{\alpha }}{\left\{{\partial }_x\left({\rho }_{\alpha }u\right)\right\}}^2〉\\ +\ell \frac{\tau }{{\theta }^2}〈c_{V_{\alpha }}{\rho }_{\alpha }{\left\{({\gamma }_{\alpha }-1)\theta {\partial }_{x}u+u{\partial }_x\theta -\frac{Q_{\alpha }}{2c_{V\alpha }{\rho }_{\alpha }}\right\}}^2〉+\frac{1}{\theta }〈Q_{\alpha }\left(1-\ell \frac{\tau ({\gamma }_{\alpha }-1)Q_{\alpha }}{4p_{\alpha }}\right)〉.\end{array}$$

(44)

It is derived from the balance equations for the mass of the components (19) and the internal energy of the mixture (43).

The sum of all terms on the right-hand side of Equation (44), except for the divergent first one, is the entropy production. Among these terms, the sum

$$\frac{\varkappa }{{\theta }^2}{\left({\partial }_x\theta \right)}^2+\frac{\nu }{\theta }{\left({\partial }_{x}u\right)}^2+\frac{1}{K{d}_0\theta }\langle d_{\alpha }^2\rangle +\frac{\langle Q_{\alpha }\rangle }{\theta }$$

is the Navier–Stokes–Fourier entropy production, and all its four terms are non-negative.

The remaining terms with the factor $\tau$ are relaxation ones, and the three penultimate terms on the right in Equation (44) contain quadratic terms with positive coefficients under the sign $\langle ·\rangle$ and therefore are non-negative. The last term on the right in Equation (44) is again non-negative for $\ell =0$, as well as for $\ell =1$ and under the condition $\tau ({\gamma }_{\alpha }-1)Q_{\alpha }⩽4p_{\alpha }$, $\alpha =\overline{1,K}$, or the more general condition

$$\tau 〈\frac{({\gamma }_{\alpha }-1)Q_{\alpha }^2}{4p_{\alpha }}〉⩽Q,$$

and then the total entropy production is non-negative as well.

The indicated entropy balance equation and the property of entropy production non-negativity remain valid for $\tau ⩾0$ again, where one should pass to a different form for the first relaxation term:

$$\frac{1}{\tau \theta }\langle {\rho }_{\alpha }{\widehat{w}}_{\alpha }^2\rangle =\frac{\tau }{\theta }〈\frac{1}{{\rho }_{\alpha }}{({\rho }_{\alpha }u{\partial }_{x}u+{\partial }_x{p}_{\alpha }-{\rho }_{\alpha }f)}^2〉.$$

In addition, this property remains valid in the absence of diffusion fluxes (for $d_0=0$): then, in Equation (44), one should simply omit both terms with $d_{\alpha }$.

4. A Spatial Discretization of the 1D Regularized System of Equations for the Dynamics of General Gas Mixtures

Define the uniform mesh ${\overline{\omega }}_h$ on $[-X,X]$, with the nodes $x_i=-X+ih$, $0⩽i⩽N$, and the step $h=\frac{2X}{N}$. Let ${\omega }_h={\overline{\omega }}_h\setminus \{-X,X\}$ be its internal part. Define also an auxiliary (conjugate) mesh ${\omega }_h^{*}$ with the nodes $x_{i+1/2}=(i+1/2)h$, $0⩽i⩽N-1$.

Let $H\left(\omega \right)$ be the space of functions defined on a mesh $\omega$. For functions $v\in H\left({\overline{\omega }}_h\right)$ and $y\in H\left({\omega }_h^{*}\right)$, we introduce the averages, shifts of the argument and difference quotients

$$\begin{array}{c}{\left[v\right]}_{i+1/2}=0.5(v_i+v_{i+1}),v_{-,i+1/2}=v_i,v_{+,i+1/2}=v_{i+1},\delta v_{i+1/2}=\frac{v_{i+1}-v_i}{h},\\ {\left[y\right]}_i^{*}=0.5(y_{i-1/2}+y_{i+1/2}),{\delta }^{*}{y}_i=\frac{y_{i+1/2}-y_{i-1/2}}{h}.\end{array}$$

Here, $v_i=v\left(x_i\right)$ and $y_{i+1/2}=y\left(x_{i+1/2}\right)$. Clearly, $[·],\delta$: $H\left({\overline{\omega }}_h\right)\to H\left({\omega }_h^{*}\right)$ and ${[·]}^{*},{\delta }^{*}$: $H\left({\omega }_h^{*}\right)\to H\left({\omega }_h\right)$.

Below, we need several mesh counterparts of the product rule for differentiation and formulas with averages for $v,u\in H\left({\overline{\omega }}_h\right)$ and $y\in H\left({\omega }_h^{*}\right)$ [28]

$$\begin{array}{c}\delta \left(uv\right)=\delta u·\left[v\right]+\left[u\right]\delta v,\end{array}$$

(45)

$$\begin{array}{c}{\delta }^{*}\left(y\left[v\right]\right)={\delta }^{*}y·v+{\left[y\delta v\right]}^{*},\end{array}$$

(46)

$$\begin{array}{c}{\delta }^{*}(\left[u\right]\left[v\right]-0.25h^2\delta u·\delta v)={\delta }^{*}\left[u\right]·v+u{\delta }^{*}\left[v\right],\end{array}$$

(47)

$$\begin{array}{c}{\left[y\right]}^{*}v={\left[y\left[v\right]\right]}^{*}-0.25h^2{\delta }^{*}\left(y\delta v\right),\end{array}$$

(48)

$$\begin{array}{c}\left[uv\right]=\left[u\right]\left[v\right]+0.25h^2\delta u·\delta v.\end{array}$$

(49)

Hereafter, in order to reduce the amount of parentheses, we write, for example, $\delta u·\left[v\right]=\left(\delta u\right)\left[v\right]$ (i.e., the multiplication sign · terminates the action of the previous operators on the left). All these formulas are, in fact, simple algebraic identities and can be straightforwardly verified (although some of them are nontrivial). Recall that Formula (47) follows from Formulas (46) and (48):

$${\delta }^{*}\left(\left[u\right]\left[v\right]\right)={\delta }^{*}\left[u\right]·v+{\left[\delta v·\left[u\right]\right]}^{*}={\delta }^{*}\left[u\right]·v+{\left[\delta v\right]}^{*}u+0.25h^2{\delta }^{*}(\delta v·\delta u),$$

as well as from the equality ${\left[\delta v\right]}^{*}={\delta }^{*}\left[v\right]$.

Following the single-component case [28], for the system of equations of the general gas mixture dynamics (1)–(7) in the absence of body forces (i.e., for $f_{\alpha }=0$), we construct the following semi-discrete (space-discrete) balance equations for the mass, momentum and total energy of the components

$$\begin{array}{c}{\partial }_t{\rho }_{\alpha }+{\delta }^{*}{j}_{\ell \alpha }=0,\end{array}$$

(50)

$${\partial }_t\left({\rho }_{\alpha }{u}_{\alpha }\right)+{\delta }^{*}(j_{\ell \alpha }\left[u_{\alpha }\right]+\left[p_{\alpha }\right])={\delta }^{*}{\mathsf{\Pi }}_{\ell \alpha }+S_{u,\alpha },$$

(51)

$$\begin{array}{c}{\partial }_t{E}_{\alpha }+{\delta }^{*}\left\{({\left[E_{\alpha }\right]}_2+\left[p_{\alpha }\right])(\left[u_{\alpha }\right]-w_{\ell \alpha })-0.25h^2\delta u_{\alpha }·\delta p_{\alpha }\right\}\\ ={\delta }^{*}(-q_{\ell \alpha }+{\mathsf{\Pi }}_{\ell \alpha }\left[u_{\alpha }\right])+S_{E,\alpha }+{\left[Q_{\alpha }\right]}^{*}\end{array}$$

(52)

on ${\omega }_h\times [0,T]$, for $\alpha =\overline{1,K}$. Herein, the pressure, specific internal and total energies of the components have the standard form

$$p_{\alpha }=({\gamma }_{\alpha }-1){\rho }_{\alpha }{\epsilon }_{\alpha }=R_{\alpha }{\rho }_{\alpha }{\theta }_{\alpha },{\epsilon }_{\alpha }=c_{V\alpha }{\theta }_{\alpha },E_{\alpha }=0.5{\rho }_{\alpha }{u}_{\alpha }^2+{\rho }_{\alpha }{\epsilon }_{\alpha }.$$

The following discretizations are applied for the mass flux and regularizing velocities

$$\begin{array}{c}{j}_{\ell \alpha }={\left[{\rho }_{\alpha }\right]}_{ln}(\left[u_{\alpha }\right]-w_{\ell \alpha }),\\ w_{\ell \alpha }=\ell \frac{\tau }{\left[{\rho }_{\alpha }\right]}\left[u_{\alpha }\right]\delta \left({\rho }_{\alpha }{u}_{\alpha }\right)+{\widehat{w}}_{\alpha },{\widehat{w}}_{\alpha }=\frac{\tau }{\left[{\rho }_{\alpha }\right]}(\left[{\rho }_{\alpha }\right]\left[u_{\alpha }\right]\delta u_{\alpha }+\delta p_{\alpha }),\end{array}$$

along with the viscous stress and heat flux

$$\begin{array}{c}{\mathsf{\Pi }}_{\ell \alpha }={\nu }_{\alpha }\delta u_{\alpha }+{\mathsf{\Pi }}_{\ell \alpha }^{\tau },{\mathsf{\Pi }}_{\ell \alpha }^{\tau }=\left[u_{\alpha }\right]\left[{\rho }_{\alpha }\right]{\widehat{w}}_{\alpha }+\ell \tau \left\{\left[u_{\alpha }\right]\delta p_{\alpha }+{\gamma }_{\alpha }{\left[p_{\alpha }\right]}_1\delta u_{\alpha }-({\gamma }_{\alpha }-1)Q_{\alpha }\right\},\end{array}$$

(53)

$$\begin{array}{c}{q}_{\ell \alpha }=-{\varkappa }_{\alpha }\delta {\theta }_{\alpha }+\ell q_{\alpha }^{\tau },-q_{\alpha }^{\tau }=\tau \left\{\left[{\rho }_{\alpha }\right]{\left[u_{\alpha }\right]}^2\left(\delta {\epsilon }_{\alpha }-\frac{{\left[p_{\alpha }\right]}_1}{{\left[{\rho }_{\alpha }\right]}^2}\delta {\rho }_{\alpha }\right)-Q_{\alpha }\left[u_{\alpha }\right]\right\}.\end{array}$$

(54)

The main sought functions ${\rho }_{\alpha }>0$, $u_{\alpha }$ and ${\theta }_{\alpha }>0$ as well as the functions $p_{\alpha }$, ${\epsilon }_{\alpha }$ and $E_{\alpha }$ are defined in space on the main mesh ${\overline{\omega }}_h$, whereas the functions $j_{\ell \alpha }$, $w_{\ell \alpha }$, ${\widehat{w}}_{\alpha }$, ${\mathsf{\Pi }}_{\ell \alpha }$, $q_{\ell \alpha }$, $\tau$, ${\nu }_{\alpha }$, ${\varkappa }_{\alpha }$ and $Q_{\alpha }$ are defined on the auxiliary mesh ${\omega }_h^{*}$. In addition, if, for example, $\tau =\mathcal{T}\left({\{{\rho }_{\alpha },u_{\alpha },{\theta }_{\alpha }\}}_{\alpha =1}^K\right)$, then one can set

$$\tau =\mathcal{T}\left({\{\left[{\rho }_{\alpha }\right],\left[u_{\alpha }\right],\left[{\theta }_{\alpha }\right]\}}_{\alpha =1}^K\right)\mathrm{or}\tau =\left[\mathcal{T}\left({\{{\rho }_{\alpha },u_{\alpha },{\theta }_{\alpha }\}}_{\alpha =1}^K\right)\right]\mathrm{on}{\omega }_h^{*}.$$

The above semi-discrete equations in the case $\ell =1$ represent the set of semi-discrete equations from [28] that are written for each of the mixture components and contain the additional exchange terms $S_{u,\alpha }$ and $S_{E,\alpha }$.

Notice that the first Formula (8) is most often used to discretize $w_{\ell \alpha }$ for $\ell =1$ [10,11,12], but here, its mesh counterpart does not hold, since $\delta \left({\rho }_{\alpha }{u}_{\alpha }^2\right)=\left[u_{\alpha }\right]\delta \left({\rho }_{\alpha }{u}_{\alpha }\right)+\left[{\rho }_{\alpha }{u}_{\alpha }\right]\delta u_{\alpha }$ due to Formula (45), but $\left[{\rho }_{\alpha }{u}_{\alpha }\right]\not\equiv \left[{\rho }_{\alpha }\right]\left[u_{\alpha }\right]$.

Along with the standard averages of the functions ${\rho }_{\alpha },u_{\alpha }$ and $p_{\alpha }$, this method involves non-standard averages of ${\rho }_{\alpha }$, $p_{\alpha }$, $E_{\alpha }$ and ${\epsilon }_{\alpha }$ of the form [28]

$$\begin{array}{c}{\left[{\rho }_{\alpha }\right]}_{ln}=\frac{1}{ln({\rho }_{\alpha -};{\rho }_{\alpha +})},{\left[p_{\alpha }\right]}_1=R_{\alpha }\left[{\rho }_{\alpha }\right]\left[{\theta }_{\alpha }\right],\\ {\left[E_{\alpha }\right]}_2=0.5{\left[{\rho }_{\alpha }\right]}_{ln}{u}_{\alpha -}{u}_{\alpha +}+{\left[{\rho }_{\alpha }\right]}_{ln}{\left[{\epsilon }_{\alpha }\right]}^{ln},\end{array}$$

(55)

$$\begin{array}{c}{\left[{\epsilon }_{\alpha }\right]}^{ln}=c_{V\alpha }ln\left({\frac{1}{\theta }}_{\alpha -};{\frac{1}{\theta }}_{\alpha +}\right)=c_{V\alpha }{\theta }_{\alpha -}{\theta }_{\alpha +}ln({\theta }_{\alpha -};{\theta }_{\alpha +}).\end{array}$$

(56)

Here, $ln(a;b)$ is the divided difference for the logarithmic function

$$ln(a;b)=\frac{lnb-lna}{b-a}\mathrm{for}a\ne b,ln(a;a)=\frac{1}{a},a>0,b>0.$$

(57)

The non-standard averages, like the simplest ones, are two-point and symmetric (and therefore, they have the approximation order $O\left(h^2\right)$ for twice continuously differentiable functions of the continuous argument). Note that $u_{\alpha -}{u}_{\alpha +}$ in ${\left[E_{\alpha }\right]}_2$ is an average of the geometric mean type for $u_{\alpha }^2$.

Recall that in order to avoid loss of accuracy at $b/a\approx 1$, instead of Formula (57) itself, some approximations are used in computations. For this purpose, it is convenient to apply the integral representation of $ln(a;b)$ and the corresponding numerical quadratures, in particular, the trapezoidal rule, the midpoint rule or Simpson’s rule

$$ln(a;b)={\int }_0^1\frac{1}{(1-r)a+br}dr\approx \frac{1}{2a}+\frac{1}{2b},\frac{2}{a+b},\frac{1}{6a}+\frac{4}{3(a+b)}+\frac{1}{6b}.$$

The first and second of these rules lead to the approximate equalities ${\left[{\epsilon }_{\alpha }\right]}^{ln}\approx \left[{\epsilon }_{\alpha }\right]$ and ${\left[{\rho }_{\alpha }\right]}_{ln}\approx \left[{\rho }_{\alpha }\right]$ which we apply in our numerical experiments in Section 6 below.

Note that using non-standard averages requires some caution. For example, clearly ${\left[c{\rho }_{\alpha }\right]}_{ln}=c{\left[{\rho }_{\alpha }\right]}_{ln}$ and ${\left[c{\epsilon }_{\alpha }\right]}^{ln}=c{\left[{\epsilon }_{\alpha }\right]}^{ln}$ for any $c=\mathrm{const}>0$. However, $\langle \left[{\rho }_{\alpha }\right]\rangle =\left[\rho \right]$ but $\langle {\left[{\rho }_{\alpha }\right]}_{ln}\rangle \not\equiv {\left[\rho \right]}_{ln}$ for $\rho =\langle {\rho }_{\alpha }\rangle$.

If one uses only the standard averages, then, in the semi-discrete entropy balance equation, there arise sign-indefinite mesh imbalances in the entropy production [28].

Let us first derive semi-discrete counterparts of Equations (10) and (11).

Lemma 2.

The semi-discrete balance equations for the kinetic and internal energies of the mixture components hold

$$\begin{array}{c}0.5{\partial }_t\left({\rho }_{\alpha }{u}_{\alpha }^2\right)+0.5{\delta }^{*}\left(j_{\ell \alpha }{u}_{\alpha -}{u}_{\alpha +}\right)+{\delta }^{*}\left[p_{\alpha }\right]·u_{\alpha }={\delta }^{*}{\mathsf{\Pi }}_{\ell \alpha }·u_{\alpha }+S_{u,\alpha }{u}_{\alpha },\\ {\partial }_t\left({\rho }_{\alpha }{\epsilon }_{\alpha }\right)+{\delta }^{*}\left(j_{\ell \alpha }{\left[{\epsilon }_{\alpha }\right]}^{ln}\right)=-{\delta }^{*}{q}_{\ell \alpha }-p_{\alpha }{\delta }^{*}(\left[u_{\alpha }\right]-w_{\ell \alpha })\end{array}$$

(58)

$$\begin{array}{c}+{[{\mathsf{\Pi }}_{\ell \alpha }\delta u_{\alpha }+w_{\ell \alpha }\delta p_{\alpha }+Q_{\alpha }]}^{*}+S_{E,\alpha }-S_{u,\alpha }{u}_{\alpha }\end{array}$$

(59)

on ${\omega }_h\times [0,T]$, with $\alpha =\overline{1,K}$.

Proof. 

Following [28], we multiply the semi-discrete momentum balance equation for the mixture components (51) by $u_{\alpha }$ and apply Formula (12). Using the semi-discrete mass balance equation for the components (50) and Formula (46) twice, we obtain the following equalities

$$\begin{array}{c}{\partial }_t{\rho }_{\alpha }·u_{\alpha }^2=-{\delta }^{*}{j}_{\ell \alpha }·u_{\alpha }^2=-{\delta }^{*}\left(j_{\ell \alpha }\left[u_{\alpha }^2\right]\right)+{\left[j_{\ell \alpha }\delta \left(u_{\alpha }^2\right)\right]}^{*},\end{array}$$

(60)

$$\begin{array}{c}{\delta }^{*}\left(j_{\ell \alpha }\left[u_{\alpha }\right]\right)·u_{\alpha }={\delta }^{*}\left(j_{\ell \alpha }{\left[u_{\alpha }\right]}^2\right)-{\left[j_{\ell \alpha }\left[u_{\alpha }\right]\delta u_{\alpha }\right]}^{*}.\end{array}$$

(61)

Taking into account also the elementary formulas

$${\left[u_{\alpha }\right]}^2=0.5\left[u_{\alpha }^2\right]+0.5u_{\alpha -}{u}_{\alpha +},0.5\delta \left(u_{\alpha }^2\right)=\left[u_{\alpha }\right]\delta u_{\alpha },$$

(62)

we derive the semi-discrete balance equation for the kinetic energy (58).

We subtract it from the semi-discrete balance equation for the total energy (52). By virtue of Formulas (47) and (46) applied twice, we have, respectively,

$$\begin{array}{c}{\delta }^{*}\left(\left[u_{\alpha }\right]\left[p_{\alpha }\right]-0.25h^2\delta u_{\alpha }·\delta p_{\alpha }\right)={\delta }^{*}\left[u_{\alpha }\right]·p_{\alpha }+{\delta }^{*}\left[p_{\alpha }\right]·u_{\alpha },\end{array}$$

(63)

$$\begin{array}{c}{\delta }^{*}\left(w_{\ell \alpha }\left[p_{\alpha }\right]\right)={\delta }^{*}{w}_{\ell \alpha }·p_{\alpha }+{\left[w_{\ell \alpha }\delta p_{\alpha }\right]}^{*},\end{array}$$

(64)

$$\begin{array}{c}{\delta }^{*}\left({\mathsf{\Pi }}_{\ell \alpha }\left[u_{\alpha }\right]\right)={\delta }^{*}{\mathsf{\Pi }}_{\ell \alpha }·u_{\alpha }+{\left[{\mathsf{\Pi }}_{\ell \alpha }\delta u_{\alpha }\right]}^{*}.\end{array}$$

(65)

Consequently, the semi-discrete balance equation for the internal energy (59) holds as well. □

In this lemma, the choice of averages $u_{\alpha -}{u}_{\alpha +}$ in ${\left[E_{\alpha }\right]}_2$ and the additional term $-0.25h^2\delta u_{\alpha }·\delta p_{\alpha }$ in the semi-discrete balance equation for the total energy (52) have already played their role, but the specific forms of averages ${\left[{\rho }_{\alpha }\right]}_{ln}$, ${\left[p_{\alpha }\right]}_1$ and ${\left[{\epsilon }_{\alpha }\right]}^{ln}$ have not yet been important and could be arbitrary. They will play their part in the next step of the method analysis.

The following result is a semi-discrete counterpart of the entropy balance equation for the general gas mixture (14), and it serves as the main one in this section.

Theorem 1.

For the semi-discrete method (50)–(56), the entropy balance equation for the mixture holds

$$\begin{array}{c}{\partial }_t\left(\rho s\right)+{\delta }^{*}〈j_{\ell \alpha }\left[s_{\alpha }\right]〉\\ ={\delta }^{*}〈\left[\frac{1}{{\theta }_{\alpha }}\right]\left({\varkappa }_{\alpha }\delta {\theta }_{\alpha }-\ell q_{\alpha }^{\tau })+B_{\alpha h}\right)〉+{\left[{\mathcal{P}}_h^{NS}+{\mathcal{P}}_h^{\tau }\right]}^{*}+〈\frac{1}{{\theta }_{\alpha }}(S_{E,\alpha }-S_{u,\alpha }{u}_{\alpha })〉\end{array}$$

(66)

on ${\omega }_h\times [0,T]$ (see definitions (13) of the entropies), where

$$\begin{array}{c}{B}_{\alpha h}:=R_{\alpha }{j}_{\ell \alpha }\left(1-\frac{\left[{\rho }_{\alpha }\right]}{{\left[{\rho }_{\alpha }\right]}_{ln}}\right)+c_{V\alpha }{j}_{\ell \alpha }\left(1-{\left[{\epsilon }_{\alpha }\right]}^{ln}\left[\frac{1}{{\epsilon }_{\alpha }}\right]\right)-0.25h^2\left({\mathsf{\Pi }}_{\ell \alpha }\delta u_{\alpha }+w_{\ell \alpha }\delta p_{\alpha }+Q_{\alpha }\right)\delta \frac{1}{{\theta }_{\alpha }},\\ {\mathcal{P}}_h^{NS}:=〈\frac{1}{{\theta }_{\alpha -}{\theta }_{\alpha +}}\left({\varkappa }_{\alpha }{\left(\delta {\theta }_{\alpha }\right)}^2+{\nu }_{\alpha }\left[{\theta }_{\alpha }\right]{\left(\delta u_{\alpha }\right)}^2\right)〉⩾0,\\ {\mathcal{P}}_h^{\tau }:=\langle \frac{1}{{\theta }_{\alpha -}{\theta }_{\alpha +}}\{\frac{1}{\tau }\left[{\rho }_{\alpha }\right]\left[{\theta }_{\alpha }\right]{\widehat{w}}_{\alpha }^2+\ell \tau R_{\alpha }\frac{{\left[{\theta }_{\alpha }\right]}^2}{\left[{\rho }_{\alpha }\right]}{\left(\delta \left({\rho }_{\alpha }{u}_{\alpha }\right)\right)}^2\\ +\ell \tau c_{V\alpha }\left[{\rho }_{\alpha }\right]{\left(\left[u_{\alpha }\right]\delta {\theta }_{\alpha }+({\gamma }_{\alpha }-1)\left[{\theta }_{\alpha }\right]\delta u_{\alpha }-\frac{Q_{\alpha }}{2c_{V\alpha }\left[{\rho }_{\alpha }\right]}\right)}^2+\left[{\theta }_{\alpha }\right]Q_{\alpha }\left(1-\ell \frac{\tau ({\gamma }_{\alpha }-1)Q_{\alpha }}{4{\left[p_{\alpha }\right]}_1}\right)\}\rangle .\end{array}$$

The last two terms on the right in Equation (66) compound the entropy production. The first three terms of ${\mathcal{P}}_h^{\tau }$ are non-negative, and the last term is non-negative for $\ell =0$, as well as for $\ell =1$ under the condition $\tau ({\gamma }_{\alpha }-1)Q_{\alpha }⩽4{\left[p_{\alpha }\right]}_1$, $\alpha =\overline{1,K}$. Moreover, $〈\frac{1}{{\theta }_{\alpha }}(S_{E,\alpha }-S_{u,\alpha }{u}_{\alpha })〉⩾0$ due to the last relation (9).

Proof. 

We use the reasoning from [28]. According to the formulas $\rho s=\langle {\rho }_{\alpha }{s}_{\alpha }\rangle$ and $\rho \epsilon =\langle c_{V\alpha }{\rho }_{\alpha }{\theta }_{\alpha }\rangle$, we first have

$$\begin{array}{c}{\partial }_t\left(\rho s\right)=〈{\partial }_t{\rho }_{\alpha }·s_{\alpha }+{\rho }_{\alpha }\left(-\frac{R_{\alpha }}{{\rho }_{\alpha }}{\partial }_t{\rho }_{\alpha }+\frac{c_{V\alpha }}{{\theta }_{\alpha }}{\partial }_t{\theta }_{\alpha }\right)〉\\ =〈{\partial }_t{\rho }_{\alpha }·s_{\alpha }-(R_{\alpha }+c_{V\alpha }){\partial }_t{\rho }_{\alpha }+{\partial }_t\left({\rho }_{\alpha }{\epsilon }_{\alpha }\right)·\frac{1}{{\theta }_{\alpha }}〉.\end{array}$$

Due to the semi-discrete mass balance equation for the component (50) and the expression

$$\begin{array}{c}{\delta }^{*}\left(j_{\ell \alpha }\left[s_{\alpha }\right]\right)={\delta }^{*}{j}_{\ell \alpha }·s_{\alpha }+{\left[j_{\ell \alpha }\delta s_{\alpha }\right]}^{*},\end{array}$$

(67)

see Formula (46), further we obtain

$${\partial }_t\left(\rho s\right)+{\delta }^{*}\langle j_{\ell \alpha }\left[s_{\alpha }\right]\rangle =〈{\left[j_{\ell \alpha }\delta s_{\alpha }\right]}^{*}+(R_{\alpha }+c_{V\alpha }){\delta }^{*}{j}_{\ell \alpha }+{\partial }_t\left({\rho }_{\alpha }{\epsilon }_{\alpha }\right)·\frac{1}{{\theta }_{\alpha }}〉.$$

By virtue of the formula

$$\begin{array}{c}{\left[j_{\ell \alpha }{\left[{\epsilon }_{\alpha }\right]}^{ln}\delta {\frac{1}{\theta }}_{\alpha }\right]}^{*}+{\delta }^{*}\left(j_{\ell \alpha }{\left[{\epsilon }_{\alpha }\right]}^{ln}\right)·\frac{1}{{\theta }_{\alpha }}-{\delta }^{*}\left(j_{\ell \alpha }{\left[{\epsilon }_{\alpha }\right]}^{ln}\left[\frac{1}{{\theta }_{\alpha }}\right]\right)=0,\end{array}$$

(68)

see Formula (46) again, we obtain

$$\begin{array}{c}{\partial }_t\left(\rho s\right)+{\delta }^{*}\langle j_{\ell \alpha }\left[s_{\alpha }\right]\rangle =〈{\left[j_{\ell \alpha }\left(\delta s_{\alpha }+{\left[{\epsilon }_{\alpha }\right]}^{ln}\delta \frac{1}{{\theta }_{\alpha }}\right)\right]}^{*}+\left\{{\partial }_t\left({\rho }_{\alpha }{\epsilon }_{\alpha }\right)+{\delta }^{*}\left(j_{\ell \alpha }{\left[{\epsilon }_{\alpha }\right]}^{ln}\right)\right\}\frac{1}{{\theta }_{\alpha }}〉\\ +{\delta }^{*}〈j_{\ell \alpha }\left(R_{\alpha }+c_{V\alpha }-{\left[{\epsilon }_{\alpha }\right]}^{ln}\left[\frac{1}{{\theta }_{\alpha }}\right]\right)〉.\end{array}$$

(69)

By definition (13) of $s_{\alpha }$, the following equalities hold

$$\begin{array}{c}\delta s_{\alpha }=-R_{\alpha }\delta ln{\rho }_{\alpha }+c_{V\alpha }\delta ln{\theta }_{\alpha }=-R_{\alpha }ln({\rho }_{\alpha -};{\rho }_{\alpha +})\delta {\rho }_{\alpha }+c_{V\alpha }ln({\theta }_{\alpha -};{\theta }_{\alpha +})\delta {\theta }_{\alpha }.\end{array}$$

(70)

The elementary formulas

$$\delta \frac{1}{{\theta }_{\alpha }}=-\frac{\delta {\theta }_{\alpha }}{{\theta }_{\alpha -}{\theta }_{\alpha +}},\left[\frac{1}{{\theta }_{\alpha }}\right]=\frac{\left[{\theta }_{\alpha }\right]}{{\theta }_{\alpha -}{\theta }_{\alpha +}}$$

(71)

are also valid. Using the first of them in the first term on the right in Equality (69), we transform the expression under the sign ${[·]}^{*}$:

$$\begin{array}{c}{j}_{\ell \alpha }\left(\delta s_{\alpha }+{\left[{\epsilon }_{\alpha }\right]}^{ln}\delta \frac{1}{{\theta }_{\alpha }}\right)=-{\left[{\rho }_{\alpha }\right]}_{ln}(\left[u_{\alpha }\right]-w_{\ell \alpha })R_{\alpha }ln({\rho }_{\alpha -};{\rho }_{\alpha +})\delta {\rho }_{\alpha }\\ +j_{\ell \alpha }\left(c_{V\alpha }ln({\theta }_{\alpha -};{\theta }_{\alpha +})-{\left[{\epsilon }_{\alpha }\right]}^{ln}\frac{1}{{\theta }_{\alpha -}{\theta }_{\alpha +}}\right)\delta {\theta }_{\alpha }=-R_{\alpha }(\left[u_{\alpha }\right]-w_{\ell \alpha })\delta {\rho }_{\alpha }.\end{array}$$

(72)

The latter equality is valid according to the definitions of the averages ${\left[{\rho }_{\alpha }\right]}_{ln}$ and ${\left[{\epsilon }_{\alpha }\right]}^{ln}$, and it is this equality that explains their above choice (55) and (56).

We transform the second term on the right in equality (69) under the sign $\langle ·\rangle$. By virtue of the semi-discrete balance equation for the mixture internal energy (59) and the formulas

$$\begin{array}{c}-{\delta }^{*}{q}_{\ell \alpha }·\frac{1}{{\theta }_{\alpha }}={\left[q_{\ell \alpha }\delta \frac{1}{{\theta }_{\alpha }}\right]}^{*}-{\delta }^{*}\left(q_{\ell \alpha }\left[\frac{1}{{\theta }_{\alpha }}\right]\right),\end{array}$$

(73)

$$\begin{array}{c}-p_{\alpha }{\delta }^{*}(\left[u_{\alpha }\right]-w_{\ell \alpha })·\frac{1}{{\theta }_{\alpha }}={\left[R_{\alpha }(\left[u_{\alpha }\right]-w_{\ell \alpha })\delta {\rho }_{\alpha }\right]}^{*}-{\delta }^{*}\left(R_{\alpha }(\left[u_{\alpha }\right]-w_{\ell \alpha })\left[{\rho }_{\alpha }\right]\right),\end{array}$$

(74)

see Formula (46), and also by virtue of Formula (48) for the third term on the right in Equation (59), we derive

$$\begin{array}{c}\left\{{\partial }_t\left({\rho }_{\alpha }{\epsilon }_{\alpha }\right)+{\delta }^{*}\left(j_{l\alpha }{\left[{\epsilon }_{\alpha }\right]}^{ln}\right)\right\}\frac{1}{{\theta }_{\alpha }}\\ =-{\delta }^{*}\left(q_{\ell \alpha }\left[\frac{1}{{\theta }_{\alpha }}\right]+R_{\alpha }(\left[u_{\alpha }\right]-w_{\ell \alpha })\left[{\rho }_{\alpha }\right]+0.25h^2({\mathsf{\Pi }}_{\ell \alpha }\delta u_{\alpha }+w_{\ell \alpha }\delta p_{\alpha }+Q_{\alpha })\delta \frac{1}{{\theta }_{\alpha }}\right)\\ +{\left[q_{\ell \alpha }\delta \frac{1}{{\theta }_{\alpha }}+R_{\alpha }(\left[u_{\alpha }\right]-w_{\ell \alpha })\delta {\rho }_{\alpha }+({\mathsf{\Pi }}_{\alpha }\delta u_{\alpha }+w_{\ell \alpha }\delta p_{\alpha }+Q_{\alpha })\left[\frac{1}{{\theta }_{\alpha }}\right]\right]}^{*}\\ +(S_{E,\alpha }-S_{u,\alpha }{u}_{\alpha })\frac{1}{{\theta }_{\alpha }}.\end{array}$$

(75)

We apply the operation $\langle ·\rangle$ to the last equality and use Formula (72) (as a result, the first term on the right in equality (69) under the sign $\langle ·\rangle$ cancels out with ${\left[R_{\alpha }(\left[u_{\alpha }\right]-w_{\ell \alpha })\delta {\rho }_{\alpha }\right]}^{*}$). We also apply the Formulas (71) and

$$\begin{array}{c}-R_{\alpha }(\left[u_{\alpha }\right]-w_{\ell \alpha })\left[{\rho }_{\alpha }\right]+j_{\ell \alpha }\left(R_{\alpha }+c_{V\alpha }-{\left[{\epsilon }_{\alpha }\right]}^{ln}\left[\frac{1}{{\theta }_{\alpha }}\right]\right)\\ =R_{\alpha }{j}_{\ell \alpha }\left(1-\frac{\left[{\rho }_{\alpha }\right]}{{\left[{\rho }_{\alpha }\right]}_{ln}}\right)+c_{V\alpha }{j}_{\ell \alpha }\left(1-{\left[{\epsilon }_{\alpha }\right]}^{ln}\left[\frac{1}{{\epsilon }_{\alpha }}\right]\right).\end{array}$$

(76)

Then, equality (69) can be rewritten as the semi-discrete entropy balance equation for the mixture

$$\begin{array}{c}{\partial }_t\left(\rho s\right)+{\delta }^{*}\langle j_{{\ell }_{\alpha }}\left[s_{\alpha }\right]\rangle ={\delta }^{*}〈-q_{\ell \alpha }\left[\frac{1}{{\theta }_{\alpha }}\right]+B_{h\alpha }〉\\ +〈{\left[\frac{1}{{\theta }_{\alpha -}{\theta }_{\alpha +}}\left({\varkappa }_{\alpha }{\left(\delta {\theta }_{\alpha }\right)}^2+{\nu }_{\alpha }\left[{\theta }_{\alpha }\right]{\left(\delta u_{\alpha }\right)}^2+A_{\ell \alpha }\right)\right]}^{*}〉+〈\frac{1}{{\theta }_{\alpha }}(S_{E,\alpha }-S_{u,\alpha }{u}_{\alpha })〉,\end{array}$$

(77)

where $B_{h\alpha }$ has been introduced above in the statement of the theorem and

$$\begin{array}{c}{A}_{\ell \alpha }:=-\ell q_{\alpha }^{\tau }\delta {\theta }_{\alpha }+\left({\mathsf{\Pi }}_{\ell \alpha }^{\tau }\delta u_{\alpha }+w_{\ell \alpha }\delta p_{\alpha }+Q_{\alpha }\right)\left[{\theta }_{\alpha }\right].\end{array}$$

(78)

The term $A_{\ell \alpha }$ can be represented as the sum of three quadratic terms with positive coefficients and a term proportional to $Q_{\alpha }$:

$$\begin{array}{c}{A}_{\ell \alpha }=\frac{1}{\tau }\left[{\rho }_{\alpha }\right]\left[{\theta }_{\alpha }\right]{\widehat{w}}_{\alpha }^2+\ell \tau R_{\alpha }\frac{{\left[{\theta }_{\alpha }\right]}^2}{\left[{\rho }_{\alpha }\right]}{\left(\delta \left({\rho }_{\alpha }{u}_{\alpha }\right)\right)}^2\\ +\ell \tau c_{V\alpha }\left[{\rho }_{\alpha }\right]{\left(\left[u_{\alpha }\right]\delta {\theta }_{\alpha }+({\gamma }_{\alpha }-1)\left[{\theta }_{\alpha }\right]\delta u_{\alpha }-\frac{Q_{\alpha }}{2c_{V\alpha }\left[{\rho }_{\alpha }\right]}\right)}^2+\left[{\theta }_{\alpha }\right]Q_{\alpha }\left(1-\ell \frac{\tau ({\gamma }_{\alpha }-1)Q_{\alpha }}{4{\left[p_{\alpha }\right]}_1}\right)\end{array}$$

(79)

according to [28] (where the case $\ell =1$ was considered, while the case $\ell =0$ is radically simpler). For the sake of completeness and the reader convenience, the proof is included below in Appendix A. Therefore, the formula holds

$$〈{\left[\frac{1}{{\theta }_{\alpha -}{\theta }_{\alpha +}}{A}_{\ell \alpha }\right]}^{*}〉={\left[{\mathcal{P}}_h^{\tau }\right]}^{*}$$

with ${\mathcal{P}}_h^{\tau }$ introduced above in the statement of the theorem, which completes the proof. □

Remark 1.

Theorem 1 can also be derived from the semi-discrete entropy balance equation for the individual mixture components

$$\begin{array}{c}{\partial }_t\left({\rho }_{\alpha }{s}_{\alpha }\right)+{\delta }^{*}\left(j_{{\ell }_{\alpha }}\left[s_{\alpha }\right]\right)={\delta }^{*}\left(-q_{\ell \alpha }\left[\frac{1}{{\theta }_{\alpha }}\right]+B_{h\alpha }\right)\\ +{\left[\frac{1}{{\theta }_{\alpha -}{\theta }_{\alpha +}}\left({\varkappa }_{\alpha }{\left(\delta {\theta }_{\alpha }\right)}^2+{\nu }_{\alpha }\left[{\theta }_{\alpha }\right]{\left(\delta u_{\alpha }\right)}^2+A_{\ell \alpha }\right)\right]}^{*}+\frac{1}{{\theta }_{\alpha }}(S_{E,\alpha }-S_{u,\alpha }{u}_{\alpha }),\end{array}$$

with $\alpha =\overline{1,K}$, see [28] in the case $S_{u,\alpha }=S_{E,\alpha }=0$. It is proved with the help of transformations quite similar to those performed above, using $s_{\alpha }$ instead of s. Applying the operation $\langle ·\rangle$ to it, we obtain the semi-discrete entropy balance equation for the mixture (77).

Note that the term ${\delta }^{*}\langle B_{\alpha h}\rangle$ in Equation (66) is a divergent difference imbalance. Of course, the appearance of the first two summands in $B_{\alpha h}$ is associated with the use of non-standard averages ${\left[{\rho }_{\alpha }\right]}_{ln}$ and ${\left[{\epsilon }_{\alpha }\right]}^{ln}$ instead of $\left[{\rho }_{\alpha }\right]$ and ${\left[\frac{1}{{\epsilon }_{\alpha }}\right]}^{-1}$ (the second of which, however, also does not coincide with the simplest $\left[{\epsilon }_{\alpha }\right]$). Moreover, $B_{\alpha h}=O\left(h^2\right)$ on continuous functions ${\rho }_{\alpha }>0$, $u_{\alpha }$ and ${\theta }_{\alpha }>0$ having bounded second derivatives ${\partial }_x^2{\rho }_{\alpha }$, ${\partial }_x^{2}u$ and ${\partial }_x^2{\theta }_{\alpha }$ and for bounded $\tau ,{\nu }_{\alpha }$ and $Q_{\alpha }$, $\alpha =\overline{1,K}$.

Similarly to the differential case, the derived semi-discrete entropy balance equation remains valid for $\tau ⩾0$, where one should pass to a different form for the first relaxation term in ${\mathcal{P}}_h^{\tau }$:

$$\frac{1}{\tau }\left[{\rho }_{\alpha }\right]\left[{\theta }_{\alpha }\right]{\widehat{w}}_{\alpha }^2=\tau \frac{\left[{\theta }_{\alpha }\right]}{\left[{\rho }_{\alpha }\right]}{(\left[{\rho }_{\alpha }\right]\left[u_{\alpha }\right]\delta u_{\alpha }+\delta p_{\alpha })}^2.$$

Applying the operation $\langle ·\rangle$ to the semi-discrete balance equations for the momentum and total energy of the components (51)–(52), as their corollaries, we obtain the semi-discrete balance equations for the momentum and total energy of the general mixture

$$\begin{array}{c}{\partial }_t\langle {\rho }_{\alpha }{u}_{\alpha }\rangle +{\delta }^{*}\left(\langle j_{\ell \alpha }\left[u_{\alpha }\right]\rangle +\left[p\right]\right)={\delta }^{*}{\mathsf{\Pi }}_{\ell },\end{array}$$

(80)

$$\begin{array}{c}{\partial }_{t}E+{\delta }^{*}\langle ({\left[E_{\alpha }\right]}_2+\left[p_{\alpha }\right])(\left[u_{\alpha }\right]-w_{\ell \alpha })-0.25h^2\delta u_{\alpha }\delta p_{\alpha }\rangle ={\delta }^{*}(-q_{\ell }+\langle {\mathsf{\Pi }}_{\ell \alpha }\left[u_{\alpha }\right]\rangle )+{\left[Q\right]}^{*}\end{array}$$

(81)

on ${\omega }_h\times [0,T]$, compare with the corresponding differential Equations (16) and (17). They involve the total viscous stress, pressure, heat flux and intensity of heat sources

$$\begin{array}{c}{\mathsf{\Pi }}_{\ell }=\langle {\nu }_{\alpha }\delta u_{\alpha }\rangle +{\mathsf{\Pi }}_{\ell }^{\tau },{\mathsf{\Pi }}_{\ell }^{\tau }=\langle \left[u_{\alpha }\right]\left[{\rho }_{\alpha }\right]{\widehat{w}}_{\alpha }\rangle +\ell \tau \left\{〈(\left[u_{\alpha }\right]\delta p_{\alpha }+{\gamma }_{\alpha }{\left[p_{\alpha }\right]}_1)\delta u_{\alpha }-{\gamma }_{\alpha }{Q}_{\alpha }〉+Q\right\},\\ p=\langle p_{\alpha }\rangle ,q_{\ell }=\langle q_{\ell \alpha }\rangle ,Q=\langle Q_{\alpha }\rangle .\end{array}$$

These equations are exploited in the next section.

5. A Spatial Discretization of the 1D Regularized System of Equations for the Dynamics of One-Velocity Gas Mixtures in the Presence of Diffusion Fluxes

In this section, we perform a spatial discretization of the system of equations for the one-velocity mixture dynamics (19)–(28) and (32)–(34). It will not be accomplished directly, but rather will be based on the aggregation of the above semi-discrete equations of the general mixture dynamics. Namely this non-standard approach will ultimately ensure the fulfillment of the correct semi-discrete entropy balance equation for the one-velocity mixture.

First, for simplicity, let the body force be absent (i.e., $f=0$). Relying on the semi-discrete balance equations for the mass of the components (50) and the momentum and total energy of the general mixture (80)–(81) and setting there $u_{\alpha }=u$ and ${\theta }_{\alpha }=\theta$ for $\alpha =\overline{1,K}$, we construct the following semi-discrete balance equations for the mass of the components and the momentum and total energy of the one-velocity mixture

$$\begin{array}{c}{\partial }_t{\rho }_{\alpha }+{\delta }^{*}(j_{\ell \alpha }+d_{\alpha })=0,\alpha =\overline{1,K},\end{array}$$

(82)

$$\begin{array}{c}{\partial }_t\left(\rho u\right)+{\delta }^{*}(j_{\ell }\left[u\right]+\left[p\right])={\delta }^{*}{\mathsf{\Pi }}_{\ell },\end{array}$$

(83)

$$\begin{array}{c}{\partial }_{t}E+{\delta }^{*}\left\{〈({\left[E_{\alpha }\right]}_2+\left[p_{\alpha }\right])(\left[u\right]-w_{\ell \alpha })〉-0.25h^2\delta u·\delta p\right\}={\delta }^{*}(-q_{\ell }+{\mathsf{\Pi }}_{\ell }\left[u\right])+{\left[Q\right]}^{*}\end{array}$$

(84)

on ${\omega }_h\times [0,T]$, where, additionally, there are discrete diffusion fluxes $d_{\alpha }$ between the component $\alpha$ and the remaining components of the mixture. The main sought functions ${\rho }_{\alpha }>0$, $\alpha =\overline{1,K}$, u and $\theta >0$ together with the functions $p_{\alpha }$, ${\epsilon }_{\alpha }$ and $E_{\alpha }$ are defined in space on the main mesh ${\overline{\omega }}_h$. In the equations, the pressure, specific internal and total energies of the components have the familiar form

$$p_{\alpha }=({\gamma }_{\alpha }-1){\rho }_{\alpha }{\epsilon }_{\alpha }=R_{\alpha }{\rho }_{\alpha }\theta ,{\epsilon }_{\alpha }=c_{V\alpha }\theta ,E_{\alpha }=0.5{\rho }_{\alpha }{u}^2+{\rho }_{\alpha }{\epsilon }_{\alpha },$$

(85)

together with the total density and corresponding quantities for the one-velocity mixture

$$\begin{array}{c}\rho =\langle {\rho }_{\alpha }\rangle ,p=\langle p_{\alpha }\rangle =R\rho \theta ,\epsilon =〈\frac{{\rho }_{\alpha }}{\rho }{\epsilon }_{\alpha }〉=c_V\theta ,E=\langle E_{\alpha }\rangle =0.5\rho u^2+\rho \epsilon ,\end{array}$$

(86)

see also the functions R and $c_V$ in (24).

There arise the following discretizations of the mass flux and regularizing velocities

$$\begin{array}{c}{j}_{\ell \alpha }={\left[{\rho }_{\alpha }\right]}_{ln}(\left[u\right]-w_{\ell \alpha }),j_{\ell }=\langle j_{\ell \alpha }\rangle ,\end{array}$$

(87)

$$\begin{array}{c}{w}_{\ell \alpha }=\ell \frac{\tau }{\left[{\rho }_{\alpha }\right]}\left[u\right]\delta \left({\rho }_{\alpha }u\right)+{\widehat{w}}_{\alpha },{\widehat{w}}_{\alpha }=\frac{\tau }{\left[{\rho }_{\alpha }\right]}(\left[{\rho }_{\alpha }\right]\left[u\right]\delta u+\delta p_{\alpha }),\end{array}$$

(88)

as well as of the viscous stress and the heat flux

$$\begin{array}{c}{\mathsf{\Pi }}_{\ell }=\nu \delta u+{\mathsf{\Pi }}_{\ell }^{\tau },{\mathsf{\Pi }}_{\ell }^{\tau }=\left[u\right]\langle \left[{\rho }_{\alpha }\right]{\widehat{w}}_{\alpha }\rangle +\ell \tau \left\{\left[u\right]\delta p+\langle {\gamma }_{\alpha }{\left[p_{\alpha }\right]}_1\rangle \delta u-\langle {\gamma }_{\alpha }{Q}_{\alpha }\rangle +Q\right\},\end{array}$$

(89)

$$\begin{array}{c}{q}_{\ell }=-\varkappa \delta \theta +q^d+\ell q^{\tau },-q^{\tau }=\tau \left\{{\left[u\right]}^2\left(\langle c_{V\alpha }\left[{\rho }_{\alpha }\right]\rangle \delta \theta -\left[\theta \right]\delta \langle R_{\alpha }{\rho }_{\alpha }\rangle \right)-Q\left[u\right]\right\};\end{array}$$

(90)

recall that here, $\nu =\langle {\nu }_{\alpha }\rangle$, $\varkappa =\langle {\varkappa }_{\alpha }\rangle$ and $Q=\langle Q_{\alpha }\rangle$. Clearly, the following formulas hold

$$\langle \left[{\rho }_{\alpha }\right]{\widehat{w}}_{\alpha }\rangle =\tau (\left[\rho \right]\left[u\right]\delta u+\delta p),\langle c_{V\alpha }\left[{\rho }_{\alpha }\right]\rangle =\left[c_V\rho \right],\langle R_{\alpha }{\rho }_{\alpha }\rangle =R\rho .$$

The new term $q^d$ in the heat flux appears due to taking into account the diffusion fluxes between the mixture components.

The functions $j_{\ell \alpha }$, $w_{\ell \alpha }$, ${\widehat{w}}_{\alpha }$, ${\mathsf{\Pi }}_{\ell }$, $q_{\ell }$, $\tau$, $\nu$, $\varkappa$ and $Q_{\alpha }$ are defined in space on the auxiliary mesh ${\omega }_h^{*}$. Moreover, we have

$$\begin{array}{c}{\left[p_{\alpha }\right]}_1=R_{\alpha }\left[{\rho }_{\alpha }\right]\left[\theta \right],\langle {\gamma }_{\alpha }{\left[p_{\alpha }\right]}_1\rangle =\langle {\gamma }_{\alpha }{R}_{\alpha }\left[{\rho }_{\alpha }\right]\rangle \left[\theta \right],\\ {\left[E_{\alpha }\right]}_2=0.5{\left[{\rho }_{\alpha }\right]}_{ln}{u}_{-}{u}_{+}+{\left[{\rho }_{\alpha }\right]}_{ln}{\left[{\epsilon }_{\alpha }\right]}^{ln},{\left[{\epsilon }_{\alpha }\right]}^{ln}=c_{V\alpha }\frac{{\theta }_{-}{\theta }_{+}}{{\left[\theta \right]}_{ln}}.\end{array}$$

We emphasize that, at first glance and quite standardly, the direct discretization of the simplified balance equations for the momentum and total energy (30) and (31) seems to be the most natural. However, it leads to a different method, since even ${\left[\rho \right]}_{ln}\left[u\right]\not\equiv \langle {\left[{\rho }_{\alpha }\right]}_{ln}\left[u\right]\rangle$, and, for such a method, the construction of an entropy balance equation with a non-negative entropy production is problematic.

In addition, we set

$$\begin{array}{c}-d_{\alpha }=d_0\left(\delta (K{G}_{\alpha }-G)+b_{\alpha }\delta \theta \right),\alpha =\overline{1,K},q^d=〈\left(\left[G_{\alpha }\right]+K^{-1}{b}_{\alpha }\left[\theta \right]\right)d_{\alpha }〉,\end{array}$$

(91)

where the Gibbs potentials of the components $G_{\alpha }$ are given by Formulas (34) and are defined in space on ${\overline{\omega }}_h$ and still $G=\langle G_{\alpha }\rangle$, whereas the functions $d_{\alpha },d_0,b_{\alpha }$ and $q^d$ are defined in space on ${\omega }_h^{*}$. Moreover, if $b_{\alpha }=B_{\alpha }({\rho }_1,\dots ,{\rho }_K,\theta )$, then to preserve the property $\langle b_{\alpha }\rangle =0$ on ${\omega }_h^{*}$, one should put, for example, $b_{\alpha }=\left[B_{\alpha }({\rho }_1,\dots ,{\rho }_K,\theta )\right]$ (but not, say, $b_{\alpha }=B_{\alpha }(\left[{\rho }_1\right],\dots ,\left[{\rho }_K\right],\left[\theta \right])$) on ${\omega }_h^{*}$.

Now, we present the semi-discrete counterpart of Proposition 1.

Proposition 2.

Let the functions $\rho >0$, u and $\theta >0$ be a solution to the following semi-discrete method for the regularized system of equations for a single-component gas dynamics [28]

$$\begin{array}{c}{\partial }_t\rho +{\delta }^{*}{j}_{\ell }=0,\end{array}$$

(92)

$$\begin{array}{c}{\partial }_t\left(\rho u\right)+{\delta }^{*}(j_{\ell }\left[u\right]+\left[p\right])={\delta }^{*}{\mathsf{\Pi }}_{\ell },\end{array}$$

(93)

$$\begin{array}{c}{\partial }_{t}E+{\delta }^{*}\left\{({\left[E\right]}_2+\left[p\right])(\left[u\right]-w_{\ell })-0.25h^2\delta u·\delta p\right\}={\delta }^{*}(-q_{\ell }+{\mathsf{\Pi }}_{\ell }\left[u\right])+{\left[Q\right]}^{*}\end{array}$$

(94)

on ${\omega }_h\times [0,T]$, where $p=(\gamma -1)\rho \epsilon =R\rho \theta$, $\epsilon =c_V\theta$ and $E=0.5\rho u^2+\rho \epsilon$ with constant $R>0$ and $c_V>0$ together with

$$\begin{array}{c}{j}_{\ell }={\left[\rho \right]}_{ln}(\left[u\right]-w_{\ell }),\end{array}$$

(95)

$$\begin{array}{c}{w}_{\ell }=\ell \frac{\tau }{\left[\rho \right]}\left[u\right]\delta \left(\rho u\right)+\widehat{w},\widehat{w}=\frac{\tau }{\left[\rho \right]}(\left[\rho \right]\left[u\right]\delta u+\delta p),\end{array}$$

(96)

$$\begin{array}{c}{\mathsf{\Pi }}_{\ell }=\nu \delta u+{\mathsf{\Pi }}_{\ell }^{\tau },{\mathsf{\Pi }}_{\ell }^{\tau }=\left[u\right]\left[\rho \right]\widehat{w}+\ell \tau \{\left[u\right]\delta p+\gamma {\left[p\right]}_1\delta u-(\gamma -1)Q\},\end{array}$$

(97)

$$\begin{array}{c}{q}_{\ell }=-\varkappa \delta \theta +\ell \left\{{\left[u\right]}^2\left(c_V\left[\rho \right]\delta \theta -R\left[\theta \right]\delta \rho \right)-Q\left[u\right]\right\},\end{array}$$

(98)

with, for example, $\tau =\left[\mathcal{T}\right(\rho ,\epsilon ,u\left)\right]>0$, $\nu =\left[\mathcal{N}\right(\rho ,\epsilon ,u\left)\right]⩾0$ and $\varkappa =\left[\mathcal{K}\right(\rho ,\epsilon ,u\left)\right]⩾0$ in space on ${\omega }_h^{*}$.

Then, for any constant $0<C_{\alpha }<1$, $\alpha =\overline{1,K}$, such that $\langle C_{\alpha }\rangle =1$, the functions ${\rho }_{\alpha }=C_{\alpha }\rho$, $\alpha =\overline{1,K}$, u and $\theta >0$ are a solution to the semi-discrete method for the regularized system of equations for the one-velocity mixture dynamics (82)–(90) with $d_{\alpha }=0$, ${\gamma }_{\alpha }=\gamma$, $R_{\alpha }=R$ and $c_{V\alpha }=c_V$, $\alpha =\overline{1,K}$, as well as $\tau =\left[\mathcal{T}\right(\rho ,\epsilon ,u\left)\right]>0$, $\nu =\left[\mathcal{N}\right(\rho ,\epsilon ,u\left)\right]⩾0$ and $\varkappa =\left[\mathcal{K}\right(\rho ,\epsilon ,u\left)\right]⩾0$ in space on ${\omega }_h^{*}$ with ρ and ε given by Formulas (86).

Proof. 

Under the hypotheses of the proposition, Formulas (85) and (86) are valid. We also have ${\widehat{w}}_{\alpha }=\widehat{w}$ and $w_{\ell \alpha }=w_{\ell }$ in expressions (88), see definitions (96); furthermore, definitions (97) and (98) coincide with the corresponding formulas for the mixture (89) and (90). Therefore, the semi-discrete balance equations for the momentum and total energy (93) and (94) coincide with the corresponding balance equations for the mixture (83) and (84). Moreover, the multiplication of the semi-discrete mass balance Equation (92) by $C_{\alpha }$ leads to the corresponding Equation (82). □

This proposition is useful for checking various properties of solutions to the constructed semi-discrete equations of the one-velocity mixture dynamics and also for testing the codes that implement them.

Let us derive semi-discrete counterparts of Equations (41)–(43).

Lemma 3.

The following semi-discrete balance equations for the total mass and kinetic and internal energies of the one-velocity mixture hold

$$\begin{array}{c}{\partial }_t\rho +{\delta }^{*}{j}_l=0,\end{array}$$

(99)

$$\begin{array}{c}0.5{\partial }_t\left(\rho u^2\right)+0.5{\delta }^{*}\left(j_{\ell }{u}_{-}{u}_{+}\right)+{\delta }^{*}\left[p\right]·u={\delta }^{*}{\mathsf{\Pi }}_{\ell }·u,\end{array}$$

(100)

$$\begin{array}{c}{\partial }_t\left(\rho \epsilon \right)+{\delta }^{*}\langle j_{\ell \alpha }{\left[{\epsilon }_{\alpha }\right]}^{ln}\rangle =-{\delta }^{*}{q}_{\ell }-\langle p_{\alpha }{\delta }^{*}(\left[u\right]-w_{\ell \alpha })\rangle +{[{\mathsf{\Pi }}_{\ell }\delta u+\langle w_{\ell \alpha }\delta p_{\alpha }\rangle +Q]}^{*}\end{array}$$

(101)

on ${\omega }_h\times [0,T]$.

Proof. 

Applying operation $\langle ·\rangle$ to the semi-discrete mass balance equation for the components (82) leads to the semi-discrete total mass balance Equation (99), due to the property $\langle d_{\alpha }\rangle =0$.

Deriving Equations (100) and (101) is similar to the proof of Lemma 2. We multiply the semi-discrete momentum balance Equation (83) by u, and we use Formulas (12) and (60)–(62), where we omit the index $\alpha$ of ${\rho }_{\alpha },u_{\alpha }$ and $j_{\ell \alpha }$, as well as the already derived Equation (99). As a result, we derive the semi-discrete balance equation for the kinetic energy of the mixture (100).

We subtract it from the semi-discrete balance equation for the total energy (84) and obtain

$$\begin{array}{c}{\partial }_t\left(\rho \epsilon \right)+{\delta }^{*}\{\langle j_{\ell \alpha }{\left[{\epsilon }_{\alpha }\right]}^{ln}\rangle +\left[p\right]\left[u\right]-\langle \left[p_{\alpha }\right]w_{\ell \alpha }\rangle -0.25h^2\delta u·\delta p\}-{\delta }^{*}\left[p\right]·u\\ ={\delta }^{*}(-q_{\ell }+{\mathsf{\Pi }}_{\ell }\left[u\right])-{\delta }^{*}{\mathsf{\Pi }}_{\ell }·u+{\left[Q\right]}^{*}.\end{array}$$

We apply Formulas (63) and (65), where we omit the index $\alpha$ of $u_{\alpha },p_{\alpha }$ and ${\mathsf{\Pi }}_{\ell \alpha }$, as well as Formula (64), and prove the semi-discrete balance equation for the internal energy of the mixture (101). □

The next result is a semi-discrete counterpart of the balance equation for the entropy of the one-velocity gas mixture (44), and it serves as the main one in this section.

Theorem 2.

Let $d_0>0$. For the semi-discrete method (82)–(91), the balance equation for the entropy of the one-velocity mixture hold

$$\begin{array}{c}{\partial }_t\left(\rho s\right)+{\delta }^{*}〈j_{\ell \alpha }\left[s_{\alpha }\right]〉={\delta }^{*}\left((\varkappa \delta \theta -\ell q^{\tau })\left[\frac{1}{\theta }\right]-\frac{\langle b_{\alpha }{d}_{\alpha }\rangle }{K}\frac{{\left[\theta \right]}^2}{{\theta }_{-}{\theta }_{+}}+B_h^{\left(d\right)}\right)+{\left[{\mathcal{P}}_h^{NS}+{\mathcal{P}}_h^{\tau }\right]}^{*}\end{array}$$

(102)

on ${\omega }_h\times [0,T]$, where

$$\begin{array}{c}{B}_h^{\left(d\right)}:=〈R_{\alpha }{j}_{\ell \alpha }\left(1-\frac{\left[{\rho }_{\alpha }\right]}{{\left[{\rho }_{\alpha }\right]}_{ln}}\right)+c_{V\alpha }{j}_{\ell \alpha }\left(1-{\left[{\epsilon }_{\alpha }\right]}^{ln}\left[\frac{1}{{\epsilon }_{\alpha }}\right]\right)〉\\ -0.25h^2\left({\mathsf{\Pi }}_{\ell }\delta u-\langle d_{\alpha }\delta G_{\alpha }\rangle +\langle w_{\ell \alpha }\delta p_{\alpha }\rangle +Q\right)\delta \frac{1}{\theta },\\ {\mathcal{P}}_h^{NS}:=\frac{1}{{\theta }_{-}{\theta }_{+}}\left\{\varkappa {\left(\delta \theta \right)}^2+\nu \left[\theta \right]{\left(\delta u\right)}^2+\frac{\left[\theta \right]}{K{d}_0}\langle d_{\alpha }^2\rangle \right\}⩾0,\\ {\mathcal{P}}_h^{\tau }:=\frac{1}{{\theta }_{-}{\theta }_{+}}\{\frac{\left[\theta \right]}{\tau }〈\left[{\rho }_{\alpha }\right]{\widehat{w}}_{\alpha }^2〉+\ell \tau {\left[\theta \right]}^2〈\frac{R_{\alpha }}{\left[{\rho }_{\alpha }\right]}{\left(\delta \left({\rho }_{\alpha }u\right)\right)}^2〉\\ +\ell \tau 〈c_{V\alpha }\left[{\rho }_{\alpha }\right]{\left(\left[u\right]\delta \theta +({\gamma }_{\alpha }-1)\left[\theta \right]\delta u-\frac{Q_{\alpha }}{2c_{V\alpha }\left[{\rho }_{\alpha }\right]}\right)}^2〉+\left[\theta \right]〈Q_{\alpha }\left(1-\ell \frac{\tau ({\gamma }_{\alpha }-1)Q_{\alpha }}{4{\left[p_{\alpha }\right]}_1}\right)〉\}.\end{array}$$

The term ${\left[{\mathcal{P}}_h^{NS}+{\mathcal{P}}_h^{\tau }\right]}^{*}$ in Equation (102) is the semi-discrete entropy production. The first three terms of ${\mathcal{P}}_h^{\tau }$ are non-negative, and the last term is non-negative for $\ell =0$, as well as for $\ell =1$ under the condition

$$\tau 〈\frac{({\gamma }_{\alpha }-1)Q_{\alpha }^2}{4{\left[p_{\alpha }\right]}_1}〉⩽Q.$$

This condition is certainly true provided that $\tau ({\gamma }_{\alpha }-1)Q_{\alpha }⩽4{\left[p_{\alpha }\right]}_1$, $\alpha =\overline{1,K}$.

Proof. 

The derivation of Equation (102) follows the same plan as in Theorem 1 with an additional transformation of the terms containing $d_{\alpha }$ and $q^d$.

Due to formulas $\rho s=\langle {\rho }_{\alpha }{s}_{\alpha }\rangle$ and $\rho \epsilon =\langle c_{V\alpha }{\rho }_{\alpha }\rangle \theta$, we have

$${\partial }_t\left(\rho s\right)=〈{\partial }_t{\rho }_{\alpha }·s_{\alpha }+{\rho }_{\alpha }\left(-\frac{R_{\alpha }}{{\rho }_{\alpha }}{\partial }_t{\rho }_{\alpha }+\frac{c_{V\alpha }}{\theta }{\partial }_t\theta \right)〉=\langle {\partial }_t{\rho }_{\alpha }·s_{\alpha }-(R_{\alpha }+c_{V\alpha }){\partial }_t{\rho }_{\alpha }\rangle +{\partial }_t\left(\rho \epsilon \right)·\frac{1}{\theta }.$$

By virtue of the semi-discrete mass balance equation for the components (82) and Formula (67), we obtain

$${\partial }_t\left(\rho s\right)+{\delta }^{*}\langle j_{\ell \alpha }\left[s_{\alpha }\right]\rangle =\langle {\left[j_{\ell \alpha }\delta s_{\alpha }\right]}^{*}\rangle +\langle (R_{\alpha }+c_{V\alpha }){\delta }^{*}{j}_{\ell \alpha }\rangle +\langle (-s_{\alpha }+R_{\alpha }+c_{V\alpha }){\delta }^{*}{d}_{\alpha }\rangle +{\partial }_t\left(\rho \epsilon \right)·\frac{1}{\theta }.$$

By virtue of Formulas (68), where we omit the index $\alpha$ of ${\theta }_{\alpha }$, together with $-s_{\alpha }+R_{\alpha }+c_{V\alpha }=G_{\alpha }/\theta$, see (34), we can write

$$\begin{array}{c}{\partial }_t\left(\rho s\right)+{\delta }^{*}\langle j_{\ell \alpha }\left[s_{\alpha }\right]\rangle =〈{\left[j_{\ell \alpha }\left(\delta s_{\alpha }+{\left[{\epsilon }_{\alpha }\right]}^{ln}\delta \frac{1}{\theta }\right)\right]}^{*}〉+\left\{{\partial }_t\left(\rho \epsilon \right)+{\delta }^{*}\langle j_{\ell \alpha }{\left[{\epsilon }_{\alpha }\right]}^{ln}\rangle \right\}\frac{1}{\theta }\\ +{\delta }^{*}〈j_{\ell \alpha }\left(R_{\alpha }+c_{V\alpha }-{\left[{\epsilon }_{\alpha }\right]}^{ln}\left[\frac{1}{\theta }\right]\right)〉+〈\frac{G_{\alpha }}{\theta }{\delta }^{*}{d}_{\alpha }〉.\end{array}$$

(103)

With the help of Formulas (70) and (71), where we omit the index $\alpha$ of ${\theta }_{\alpha }$, in the first term on the right in equality (103), the expression under the sign ${[·]}^{*}$ can be transformed similarly to equalities (72) as

$$\begin{array}{c}{j}_{\ell \alpha }\left(\delta s_{\alpha }+{\left[{\epsilon }_{\alpha }\right]}^{ln}\delta \frac{1}{\theta }\right)=-R_{\alpha }(\left[u\right]-w_{\ell \alpha })\delta {\rho }_{\alpha }.\end{array}$$

(104)

Let us transform the second term on the right in equality (103). Due to the semi-discrete balance equation for the internal energy of the mixture (101), as well as Formulas (73) and (74), where we omit the index $\alpha$ of $q_{\ell \alpha },{\theta }_{\alpha }$ and $u_{\alpha }$, together with Formula (48), we obtain

$$\begin{array}{c}\left\{{\partial }_t\left(\rho \epsilon \right)+{\delta }^{*}\langle j_{l\alpha }{\left[{\epsilon }_{\alpha }\right]}^{ln}\rangle \right\}\frac{1}{\theta }\\ =-{\delta }^{*}\left\{q_{\ell }\left[\frac{1}{\theta }\right]+\langle R_{\alpha }(\left[u\right]-w_{\ell \alpha })\left[{\rho }_{\alpha }\right]\rangle +0.25h^2({\mathsf{\Pi }}_{\ell }\delta u+\langle w_{\ell \alpha }\delta p_{\alpha }\rangle +Q)\delta \frac{1}{\theta }\right\}\\ +\langle {\left[R_{\alpha }(\left[u\right]-w_{\ell \alpha })\delta {\rho }_{\alpha }\right]}^{*}\rangle +{\left[q_{\ell }\delta \frac{1}{\theta }+({\mathsf{\Pi }}_{\ell }\delta u+\langle w_{\ell \alpha }\delta p_{\alpha }\rangle +Q)\left[\frac{1}{\theta }\right]\right]}^{*}\end{array}$$

by analogy with equality (75).

By virtue of Formula (104), the first term on the right in equality (103) cancels out the second term on the right in the last formula. Now, due to the definitions of $q_{\ell }$ and ${\mathsf{\Pi }}_{\ell }$, see expressions (90) and (89), and using Formulas (76) and (71), where we omit the index $\alpha$ of $u_{\alpha }$ and ${\theta }_{\alpha }$, equality (103) can be rewritten as the semi-discrete balance equation for the mixture entropy

$$\begin{array}{c}{\partial }_t\left(\rho s\right)+{\delta }^{*}\langle j_{\ell \alpha }\left[s_{\alpha }\right]\rangle ={\delta }^{*}\left(\left[\frac{1}{\theta }\right](\varkappa \delta \theta -\ell q^{\tau })+B_h\right)+{\left[\frac{1}{{\theta }_{-}{\theta }_{+}}\left\{\varkappa {\left(\delta \theta \right)}^2+\nu \left[\theta \right]{\left(\delta u\right)}^2+A_{\ell }\right\}\right]}^{*}\\ -{\delta }^{*}\left(q^d\left[\frac{1}{\theta }\right]\right)+{\left[q^d\delta \frac{1}{\theta }\right]}^{*}+〈\frac{G_{\alpha }}{\theta }{\delta }^{*}{d}_{\alpha }〉,\end{array}$$

(105)

where

$$\begin{array}{c}{B}_h:=〈R_{\alpha }{j}_{\ell \alpha }\left(1-\frac{\left[{\rho }_{\alpha }\right]}{{\left[{\rho }_{\alpha }\right]}_{ln}}\right)+c_{V\alpha }{j}_{\ell \alpha }\left(1-{\left[{\epsilon }_{\alpha }\right]}^{ln}\left[\frac{1}{{\epsilon }_{\alpha }}\right]\right)〉-0.25h^2({\mathsf{\Pi }}_{\ell }\delta u+\langle w_{\ell \alpha }\delta p_{\alpha }\rangle +Q)\delta \frac{1}{\theta },\\ A_{\ell }:=-\ell q^{\tau }\delta \theta +\left({\mathsf{\Pi }}_{\ell }^{\tau }\delta u+\langle w_{\ell \alpha }\delta p_{\alpha }\rangle +Q\right)\left[\theta \right]=\langle A_{\ell \alpha }\rangle ,\\ A_{\ell \alpha }:=-\ell q_{\alpha }^{\tau }\delta \theta +\left({\mathsf{\Pi }}_{\ell \alpha }^{\tau }\delta u+w_{\ell \alpha }\delta p_{\alpha }+Q_{\alpha }\right)\left[\theta \right],\end{array}$$

with

$$\begin{array}{c}{q}_{\alpha }^{\tau }=\tau \left\{{\left[u\right]}^2\left(c_{V\alpha }\left[{\rho }_{\alpha }\right]\delta \theta -R_{\alpha }\left[\theta \right]\delta {\rho }_{\alpha }\right)-Q_{\alpha }\left[u\right]\right\},\\ {\mathsf{\Pi }}_{\ell \alpha }^{\tau }=\left[u\right]\left[{\rho }_{\alpha }\right]{\widehat{w}}_{\alpha }+\ell \tau \left\{\left[u\right]\delta p_{\alpha }+{\gamma }_{\alpha }{\left[p_{\alpha }\right]}_1\delta u-({\gamma }_{\alpha }-1)Q_{\alpha }\right\}.\end{array}$$

Let us transform the terms with $q^d$ in the balance Equation (105) with the help of Formula (46):

$$\begin{array}{c}-{\delta }^{*}\left(q^d\left[\frac{1}{\theta }\right]\right)+{\left[q^d\delta \frac{1}{\theta }\right]}^{*}+〈\frac{G_{\alpha }}{\theta }{\delta }^{*}{d}_{\alpha }〉\\ ={\delta }^{*}\left(〈d_{\alpha }\left[\frac{G_{\alpha }}{\theta }\right]〉-q^d\left[\frac{1}{\theta }\right]\right)+{\left[q^d\delta \frac{1}{\theta }-〈d_{\alpha }\delta \frac{G_{\alpha }}{\theta }〉\right]}^{*}.\end{array}$$

(106)

For the term on the right under the sign of ${\delta }^{*}$, with the help of Formula (49), the definition $q^d=\langle (\left[G_{\alpha }\right]+K^{-1}{b}_{\alpha }\left[\theta \right])d_{\alpha }\rangle$ and the second Formula (71) with $\theta$ in the role of ${\theta }_{\alpha }$, we have

$$\begin{array}{c}〈d_{\alpha }\left[\frac{G_{\alpha }}{\theta }\right]〉-q^d\left[\frac{1}{\theta }\right]=〈d_{\alpha }\left(\left[G_{\alpha }\right]\left[\frac{1}{\theta }\right]+0.25h^2\delta G_{\alpha }·\delta \frac{1}{\theta }\right)〉-\langle (\left[G_{\alpha }\right]+K^{-1}{b}_{\alpha }\left[\theta \right])d_{\alpha }\rangle \left[\frac{1}{\theta }\right]\\ =0.25h^2\langle d_{\alpha }\delta G_{\alpha }\rangle \delta \frac{1}{\theta }-\frac{\langle b_{\alpha }{d}_{\alpha }\rangle {\left[\theta \right]}^2}{K{\theta }_{-}{\theta }_{+}}.\end{array}$$

For the term on the right in (106) under the sign of ${[·]}^{*}$, with the help of formula (45), the definition of $q^d$ and formulas (71) with $\theta$ in the role of ${\theta }_{\alpha }$, as well as the definition of $d_{\alpha }$ (see expression (91)) and the property $\langle d_{\alpha }\rangle =0$, we also obtain

$$\begin{array}{c}{q}^d\delta \frac{1}{\theta }-〈d_{\alpha }\delta \frac{G_{\alpha }}{\theta }〉=\left(q^d-\langle d_{\alpha }\left[G_{\alpha }\right]\rangle \right)\delta \frac{1}{\theta }-\langle d_{\alpha }\delta G_{\alpha }\rangle \left[\frac{1}{\theta }\right]\\ =-\frac{\left[\theta \right]}{K{\theta }_{-}{\theta }_{+}}〈d_{\alpha }\left(b_{\alpha }\delta \theta +\delta (K{G}_{\alpha }-G)+\delta G\right)〉=\frac{\left[\theta \right]}{K{\theta }_{-}{\theta }_{+}}\frac{\langle d_{\alpha }^2\rangle }{d_0}.\end{array}$$

Now, the semi-discrete balance equation for the mixture entropy (105) takes the form

$$\begin{array}{c}{\partial }_t\left(\rho s\right)+{\delta }^{*}〈j_{\ell \alpha }\left[s_{\alpha }\right]〉={\delta }^{*}\left((\varkappa \delta \theta -\ell q^{\tau })\left[\frac{1}{\theta }\right]-\frac{\langle b_{\alpha }{d}_{\alpha }\rangle }{K}\frac{{\left[\theta \right]}^2}{{\theta }_{-}{\theta }_{+}}+B_h+0.25h^2\langle d_{\alpha }\delta G_{\alpha }\rangle \delta \frac{1}{\theta }\right)\\ +{\left[\frac{1}{{\theta }_{-}{\theta }_{+}}\left\{\varkappa {\left(\delta \theta \right)}^2+\nu \left[\theta \right]{\left(\delta u\right)}^2+\frac{\left[\theta \right]}{K{d}_0}\langle d_{\alpha }^2\rangle +\langle A_{\ell \alpha }\rangle \right\}\right]}^{*}.\end{array}$$

By virtue of Formula (79) with $u_{\alpha }=u$ and ${\theta }_{\alpha }=\theta$ for $\alpha =\overline{1,K}$, the term $A_{\ell \alpha }$ can be represented as the following sum

$$\begin{array}{c}{A}_{\ell \alpha }=\frac{\left[\theta \right]}{\tau }\left[{\rho }_{\alpha }\right]{\widehat{w}}_{\alpha }^2+\ell \tau {\left[\theta \right]}^2\frac{R_{\alpha }}{\left[{\rho }_{\alpha }\right]}{\left(\delta \left({\rho }_{\alpha }u\right)\right)}^2\\ +\ell \tau c_{V\alpha }\left[{\rho }_{\alpha }\right]{\left(\left[u\right]\delta \theta +({\gamma }_{\alpha }-1)\left[\theta \right]\delta u-\frac{Q_{\alpha }}{2c_{V\alpha }\left[{\rho }_{\alpha }\right]}\right)}^2+\left[\theta \right]Q_{\alpha }\left(1-\ell \frac{\tau ({\gamma }_{\alpha }-1)Q_{\alpha }}{4{\left[p_{\alpha }\right]}_1}\right),\end{array}$$

that after applying the operation $\langle ·\rangle$ allows us to complete the proof. □

As in the differential case, the indicated property of the entropy production non-negativity remains valid for $\tau ⩾0$, where one should pass to a different form for the first relaxation term in ${\mathcal{P}}_h^{\tau }$:

$$\frac{\left[\theta \right]}{\tau }〈\left[{\rho }_{\alpha }\right]{\widehat{w}}_{\alpha }^2〉=\tau \left[\theta \right]〈\frac{1}{\left[{\rho }_{\alpha }\right]}{(\left[{\rho }_{\alpha }\right]\left[u\right]\delta u+\delta p_{\alpha })}^2〉.$$

In addition, this property is also valid in the absence of diffusion fluxes (for $d_0=0$), when in the balance Equation (102), one should simply omit all three terms with $d_{\alpha }$, including such terms in $B_h^{\left(d\right)}$ and ${\mathcal{P}}_h^{NS}$.

To complete the section, following [28], we generalize the constructed semi-discrete method and the obtained results to the case of any f. Let us generalize the right-hand sides of the semi-discrete momentum and total energy balance Equations (83) and (84):

$$\begin{array}{c}{\partial }_t\left(\rho u\right)+{\delta }^{*}(j_{\ell }\left[u\right]+\left[p\right])=\delta {\mathsf{\Pi }}_{\ell }+{\left[{\rho }_{*\ell }f\right]}^{*},\\ {\partial }_{t}E+{\delta }^{*}\left\{〈({\left[E_{\alpha }\right]}_2+\left[p_{\alpha }\right])(\left[u\right]-w_{\ell \alpha })〉-0.25h^2\delta u·\delta p\right\}\\ ={\delta }^{*}(-q_{\ell }+{\mathsf{\Pi }}_{\ell }\left[u\right])+{\left[Q\right]}^{*}+〈{\left[\left[{\rho }_{\alpha }\right](\left[u\right]-w_{\ell \alpha })f\right]}^{*}〉+0.25h^2{\left[{\rho }_{*\ell }\delta u·\delta \frac{1}{\theta }·f\right]}^{*}\theta ,\end{array}$$

(107)

where the functions ${\rho }_{*\ell }:=\left[\rho \right]-\ell \tau \delta \left(\rho u\right)$ and f are defined in space on ${\omega }_h^{*}$. We also generalize expression (88) for ${\widehat{w}}_{\alpha }$ as follows

$$\begin{array}{c}{\widehat{w}}_{\alpha }=\frac{\tau }{\left[{\rho }_{\alpha }\right]}\left(\left[{\rho }_{\alpha }\right]\left[u\right]\delta u+\delta p_{\alpha }-\left[{\rho }_{\alpha }\right]f\right).\end{array}$$

(108)

On the right-hand side of the semi-discrete balance equation for the internal energy of the mixture (101), the above expressions with f generate the following additional term

$$\mathsf{\Psi }=〈{\left[\left[{\rho }_{\alpha }\right](\left[u\right]-w_{\ell \alpha })f\right]}^{*}〉-{\left[{\rho }_{*\ell }f\right]}^{*}u+0.25h^2{\left[{\rho }_{*\ell }\delta u·\delta \frac{1}{\theta }·f\right]}^{*}\theta .$$

(109)

When deriving the semi-discrete balance equation for the mixture entropy (102), one should multiply this term by $\frac{1}{\theta }$ and apply the formulas

$$\begin{array}{c}{\left[\left[{\rho }_{\alpha }\right](\left[u\right]-w_{\ell \alpha })f\right]}^{*}\frac{1}{\theta }={\left[\left[{\rho }_{\alpha }\right](\left[u\right]-w_{\ell \alpha })f\left[\frac{1}{\theta }\right]\right]}^{*}-0.25h^2{\delta }^{*}\left(\left[{\rho }_{\alpha }\right](\left[u\right]-w_{\ell \alpha })f\delta \frac{1}{\theta }\right),\\ -{\left[{\rho }_{*\ell }f\right]}^{*}\frac{u}{\theta }=-{\left[{\rho }_{*\ell }\left[u\right]f\left[\frac{1}{\theta }\right]\right]}^{*}-{\left[{\rho }_{*\ell }\left(\left[\frac{u}{\theta }\right]-\left[u\right]\left[\frac{1}{\theta }\right]\right)f\right]}^{*}+0.25h^2{\delta }^{*}\left({\rho }_{*\ell }f\delta \frac{u}{\theta }\right),\end{array}$$

see Formula (48). Further application of Formula (49) makes it possible to cancel out the last term of $\frac{\mathsf{\Psi }}{\theta }$, see expression (109), and the second term on the right in the last formula (this is why the non-standard last term on the right in Equation (107) was added). This leads to the equalities

$$\begin{array}{c}\frac{\mathsf{\Psi }}{\theta }={\left[\{\langle \left[{\rho }_{\alpha }\right](\left[u\right]-w_{\ell \alpha })\rangle -{\rho }_{*\ell }\left[u\right]\}f\left[\frac{1}{\theta }\right]\right]}^{*}+{\delta }^{*}{C}_h=-{\left[\langle \left[{\rho }_{\alpha }\right]{\widehat{w}}_{\alpha }\rangle f\left[\frac{1}{\theta }\right]\right]}^{*}+{\delta }^{*}{C}_h\end{array}$$

with

$$\begin{array}{c}{C}_h=-0.25h^2\left(\langle \left[{\rho }_{\alpha }\right](\left[u\right]-w_{\ell \alpha })\rangle \delta \frac{1}{\theta }-{\rho }_{*\ell }\delta \frac{u}{\theta }\right)f=0.25h^2\left(\langle \left[{\rho }_{\alpha }\right]{\widehat{w}}_{\alpha }\rangle \delta \frac{1}{\theta }+{\rho }_{*\ell }\delta u·\left[\frac{1}{\theta }\right]\right)f,\end{array}$$

where two formulas

$$\langle \left[{\rho }_{\alpha }\right](\left[u\right]-w_{\ell \alpha })\rangle -{\rho }_{*\ell }\left[u\right]=-\langle \left[{\rho }_{\alpha }\right]{\widehat{w}}_{\alpha }\rangle ,\delta \frac{u}{\theta }=\left[u\right]\delta \frac{1}{\theta }+\delta u·\left[\frac{1}{\theta }\right]$$

have additionally been used. As a result, the additional term $-\left[{\rho }_{\alpha }\right]{\widehat{w}}_{\alpha }f$ appears in $A_{\ell \alpha }$, and the summand $C_h$ should be added to $B_h^{\left(d\right)}$. With such f-dependent extensions (including the generalized Formula (108) for ${\widehat{w}}_{\alpha }$), Theorem 2 remains valid.

Notice that $C_h=O\left(h^2\right)$ (the same applies to the last term in Equation (107)) on continuous functions $\rho$, u and $\theta >0$ having bounded derivatives ${\partial }_x{\rho }_{\alpha }$, ${\partial }_{x}u$ and ${\partial }_x\theta$, and for bounded $\tau$ and f.

6. Numerical Experiments

We consider the binary one-velocity mixtures and three test examples well known in the literature, using the following piecewise constant initial data $({\rho }_1,{\rho }_2,p,u,\gamma ){|}_{t=0}=({\rho }_1^0,{\rho }_2^0,p^0,u^0,{\gamma }^0)$, with a discontinuity between two gases:

$$\begin{array}{c}({\rho }_1^0,{\rho }_2^0,p^0,u^0,{\gamma }^0)\left(x\right)=\left\{\begin{array}{cc}({\rho }_{1l},{\rho }_{2l},p_l,u_l,{\gamma }_l),\hfill & x<0\hfill \\ ({\rho }_{1r},{\rho }_{2r},p_r,u_r,{\gamma }_r),\hfill & x⩾0\hfill \end{array}\right.;\end{array}$$

moreover, we have ${\rho }_{1r}={\rho }_{2l}=0$ (although, instead of that, in computations, we set them very small, namely equal ${10}^{-8}$). Their parameters to the left and right of the discontinuity at $x=0$ and the final time of computations $t_{fin}$ are given in

Table 1. The parameters of the initial data to the left and right of the discontinuity between two gases and the final time of computations.

Example	$\mathit{\rho }$	p	u	$\mathit{\gamma }$	${\mathit{t}}_{\mathit{fin}}$
(1) left	1	1	0	1.4	0.2
(1) right	0.125	0.1	0	1.6
(2) left	1	2	0	2	0.2
(2) right	0.125	0.1	0	1.4
(3) left	1	500	0	1.4	0.011
(3) right	1	0.2	0	1.6

. Note that, in Examples 1–3, the initial pressure drop increases and, respectively, equals $\frac{p_l}{p_r}=10,20,2500$, which amplifies computational complexity. The initial temperature ${\theta }^0$ is calculated according to the above given formulas $p^0=R{\rho }^0{\theta }^0=(({\gamma }_l-1)c_{Vl}{\rho }_{1l}+({\gamma }_r-1)c_{Vr}{\rho }_{2r}){\theta }^0$, with the simplest choice of $c_{Vl}=c_{Vr}=1$.

We also set $X=0.5$ in all the examples, and the boundary values of the sought functions at $x=\pm X$ in time are set constant and the same as their values given at $t=0$.

We define a non-uniform mesh in time $0=t_0<t_1<\dots <t_{\overline{m}}=t_{fin}$, with the steps $h_{tm}=t_m-t_{m-1}$. We take the relaxation parameter and the artificial viscosity and heat conductivity coefficients in the form

$${\tau }_i^m=\frac{ah}{c_{si}^m},\mu =\tau p=\tau (p_1+p_2),\varkappa =\tau (a_{Pr,l}{\gamma }_l{c}_{Vl}{p}_1+a_{Pr,r}{\gamma }_r{c}_{Vr}{p}_2)$$

standard for the single-component case [11,12]. Here, $0<a<1$ is a parameter; recall also that $h=\frac{2X}{N}$ is the step of the spatial mesh, $c_s:=\sqrt{(\gamma -1)\gamma \epsilon }$ is the speed of sound in the mixture and, for example, $c_{si}^m=c(x_i,t_m)$. In addition, $a_{Pr,l}$ and $a_{Pr,r}$ are the inverse Prandtl numbers (which can be considered as adjusting parameters as well); we take them equal 1 except for the last computation.

We apply the discretization in space constructed in the previous Section 5 combined with the simplest explicit Euler method for the discretization in time. We use the automatic non-uniform mesh in time, with the steps $h_{tm}=t_m-t_{m-1}$ such that

$$h_{tm}=\frac{\beta h}{{max}_{0⩽i⩽N}(c_{si}^{m-1}+|u_i^{m-1}\left|\right)},1⩽m⩽\overline{m}-1,h_{t\overline{m}}=t_{fin}-t_{\overline{m}-1}⩽\frac{\beta h}{{max}_{0⩽i⩽N}\left(c_{si}^{\overline{m}-1}+\left|u_i^{\overline{m}-1}\right|\right)},$$

where $\beta$ is a parameter (the Courant-type number). Note that if there appear values of ${\rho }_1$ or ${\rho }_2$ less than ${10}^{-10}$ at the upper time level, we replace them by ${10}^{-10}$. We adjust the parameters a and $\beta$ in each example. Proposition 2 (more precisely, its natural fully discrete counterpart) was used to test the code initially.

Example 1 (moderate two-gas shock-tube problem) from ([34], Test 5.3), see also ([35], Test 3.3). The final solution at $t=t_{fin}$ has jumps in the values of ${\rho }_1$, ${\rho }_2$, $\rho$ and $\theta$ but not p and u at the contact discontinuity between the two gases, as well as a rarefaction wave in gas 1 to the left and a shock wave (the strong discontinuity, with jumps in the values of p and u as well) in gas 2 to the right of the contact discontinuity. The functions ${\rho }_1$, $\rho$ and p are non-increasing, whereas ${\rho }_2$, u and $\theta$ are non-monotone, with the maximal values of ${\rho }_2$ and u in front of the shock; also, ${\rho }_2$ is piecewise constant. The final maximal Mach number is $M_{max}:={max}_{0⩽i⩽N}\frac{|u_i^{\overline{m}}|}{c_{si}^{\overline{m}}}\approx 0.91$, so the flow is subsonic.

In the QGD case, i.e., for $\ell =1$, the results for $a=0.25$, $\beta =0.4$ and $N=1601$ are shown in Figure 1 containing the graphs of ${\rho }_1$, ${\rho }_2$, $\rho$, p, u and $\theta$. For much smaller $N=401$, they are similar; however, as usual, the graph slopes near the points of strong and especially contact discontinuities are smaller, whereas the slopes for $N=801$ are already much closer to those for $N=1601$. The results correspond well to those given in [34,35].

Moreover, in the simpler QHD case, i.e., for $\ell =0$, the results for $a=0.6$ and the same $\beta$ and N are very close, only with a very small ledge in the graph of u at the contact discontinuity and the graph of ${\rho }_2$ just to the right of it, see the same figure. Hereafter, more minor differences become visible after magnifying our figures several times.

Example 2 (a modified two-gas Sod problem) from ([36], Test 3.2). In this example, ${\gamma }_l>{\gamma }_r$ in contrast to Example 1, and also the difference $|{\gamma }_l-{\gamma }_r|$ is larger. The final solution is in general similar to Example 1, but the maximal value of $\theta$ is now in front of the shock, not at the left boundary. The final maximal Mach number is $M_{max}\approx 0.94$, so the flow is subsonic once again.

In the QGD case, the results for $a=0.25$, $\beta =0.3$ and $N=2001$ are presented in Figure 2 (for $N=501$, they are similar, but the graph slopes near the points of contact and strong discontinuities are smaller once again). They correspond well to those from [36] except for a very small hollow at the point of contact discontinuity.

In the QHD case, the results for $a=0.6$, $\beta =0.3$ and $N=2001$ are very close once again; a unique visible difference is slightly worse behavior of the graph of ${\rho }_2$ just to the right of the contact discontinuity in the same figure.

Example 3 (stiff two-gas shock-tube problem) from ([34], Test 5.4). In the original paper, the initial left and right pressures were confused, and $t_{fin}$ was not specified, so we adjusted it according to the given graphs; note that it is much less than in the previous examples. Recall that, in this example, the initial pressure drop equals 2500 and is much larger than previously. Concerning the final values, now, the support of the maximal value of ${\rho }_2$ is more narrow, the maximal value of $\rho$ is in front of the shock and $\theta$ becomes non-increasing, with large jumps in their values. In addition, $M_{max}\approx 1.44$, so the flow is partially supersonic now (recall that $u=0$ and thus $M=0$ near the boundaries).

In the QGD case, the results for $a=0.25$, $\beta =0.1$ and $N=4001$ are given in Figure 3, and they are in accordance with those in [34]. For smaller $N=2001$, the quality of graphs of ${\rho }_2$ and $\rho$ near the contact discontinuity becomes worse.

In contrast, in the QHD case, now, the graphs of ${\rho }_2$, u and p are satisfactory even for $N=1001$, but the quality of the graphs of $\rho$ and especially of ${\rho }_1$ and $\theta$ is much worse near the contact discontinuity (even for larger $a=0.6$), and we omit them. However, the situation can be radically improved by taking larger $a_{Pr,l}=a_{Pr,r}=10$ and applying another known formula for $\tau$:

$${\tau }_i^m=\frac{ah}{c_{si}^m+\left|u_i^m\right|},$$

that is also often used for transonic or supersonic flows. The results become very close and are very slightly worse only for ${\rho }_1$ and u near the contact discontinuity; see the same figure.