Compressible multi-component flow in porous media with Maxwell-Stefan diffusion

TL;DR

This paper develops a Darcy-scale model for compressible multi-component flow in porous media using Maxwell-Stefan diffusion, ensuring thermodynamic consistency and analyzing solution behavior over time.

Contribution

It introduces a new thermodynamically consistent Darcy-scale model incorporating Maxwell-Stefan diffusion for multi-component flow in porous media.

Findings

Existence of global classical solutions near equilibrium.

Solutions tend to a parabolic limit system over time.

Model captures cross-diffusive effects accurately.

Abstract

We introduce a Darcy-scale model to describe compressible multi-component flow in a fully saturated porous medium. In order to capture cross-diffusive effects between the different species correctly, we make use of the Maxwell--Stefan theory in a thermodynamically consistent way. For inviscid flow the model turns out to be a nonlinear system of hyperbolic balance laws. We show that the dissipative structure of the Maxwell-Stefan operator permits to guarantee the existence of global classical solutions for initial data close to equilibria. Furthermore it is proven by relative entropy techniques that solutions of the Darcy-scale model tend in a certain long-time regime to solutions of a parabolic limit system.

Equations

$$\begin{array}{rl}{\partial }_t\rho +\mathrm{div}(\mathbf{m})& =0,\\ {\partial }_t\mathbf{m}+\mathrm{div}\left(\frac{\mathbf{m}{\mathbf{m}}^{\top }}{\rho }+p(\rho ){\mathcal{I}}_d\right)& =-M\mathbf{m}\end{array}$$

$${\partial }_t\eta (\rho ,\mathbf{m})+\mathrm{div}\mathbf{q}(\rho ,\mathbf{m})\le 0\text{in }{\mathcal{D}}^{\prime }({\mathbb{R}}^d\times (0,T))$$

$$\eta (\rho ,\mathbf{m})=\frac{1}{2}\frac{|\mathbf{m}{|}^2}{\rho }+\rho \psi (\rho ),\mathbf{q}(\rho ,\mathbf{m})=\frac{1}{2}\mathbf{m}\frac{|\mathbf{m}{|}^2}{{\rho }^2}+\mathbf{m}(\psi (\rho )+\rho {\psi }^{\prime }(\rho )).$$

$$\rho \psi (\rho )=(\rho \psi {)}^{\prime }(\rho )-p(\rho ).$$

$$\begin{array}[]{rcl}\varepsilon\partial_{t}\rho^{\varepsilon}+\operatorname{div}(\mathbf{m}^{\varepsilon})&=&0,\\ \varepsilon\partial_{t}(\mathbf{m}^{\varepsilon})+\operatorname{div}\left(\dfrac{\mathbf{m}^{\varepsilon}\mathbf{m}^{\varepsilon\top}}{\rho^{\varepsilon}}+p(\rho^{\varepsilon})\mathcal{I}_{d}\right)&=&-\frac{1}{\varepsilon}M\mathbf{m}^{\varepsilon}\end{array}\text{ in }{\mathbb{R}}^{d}\times(0,T).$$

$${\partial }_t\bar{\rho }-M^{-1}\mathrm{div}(\nabla p(\bar{\rho }))=0\text{ in }{\mathbb{R}}^d\times (0,T).$$

$$\rho :=\sum _{i=1}^n{\rho }_i,\mathbf{v}:=\frac{1}{\rho }\sum _{i=1}^n{\rho }_i{\mathbf{v}}_i.$$

$${\mathbf{u}}_i:={\mathbf{v}}_i-\mathbf{v}\in {\mathbb{R}}^d.$$

$${\partial }_t{\rho }_i+\mathrm{div}({\rho }_i{\mathbf{v}}_i)$$

$${\partial }_t({\rho }_i{\mathbf{v}}_i)+\mathrm{div}({\rho }_i{\mathbf{v}}_i{\mathbf{v}}_i^{\top }+p_i({\rho }_i){\mathcal{I}}_d)$$

$$\sum _{i=1}^n{\mathbf{f}}_i=0.$$

$$h_i({\rho }_i)+p_i({\rho }_i)={\rho }_i{h}_i^{\prime }({\rho }_i).$$

$$p_i^{\prime }({\rho }_i)>0.$$

$$h({\rho }_1,\dots ,{\rho }_n):=\sum _{i=1}^n{h}_i({\rho }_i)$$

$$\zeta =-\sum _{i=1}^n{\mathbf{u}}_i\cdot {\mathbf{f}}_i.$$

$$-\sum _{i=1}^{n-1}({\mathbf{u}}_i-{\mathbf{u}}_n)\cdot {\mathbf{f}}_i\ge 0.$$

$$\tilde{\mathcal{T}}:=({\tau }_{ij}{)}_{i,j=1}^{n-1}\succ 0,\text{ with }{\tau }_{ij}={\tau }_{ij}({\rho }_i,{\rho }_j)$$

$${\mathbf{f}}_i=-\sum _{j=1}^{n-1}{\tau }_{ij}({\mathbf{u}}_j-{\mathbf{u}}_n),i=1,...,n-1.$$

$${\tau }_{nj}=-\sum _{i=1}^{n-1}{\tau }_{ij},j=1,...,n-1,{\tau }_{in}=-\sum _{j=1}^{n-1}{\tau }_{ij},i=1,...,n.$$

$${\tau }_{ij}\le 0\text{ for all }i\ne j.$$

$${\tau }_{ij}=-{\lambda }_{ij}({\rho }_i,{\rho }_j){\rho }_i{\rho }_j(i\ne j).$$

$$\Lambda=\Lambda(\mathbf{r})={\Big{(}\lambda_{ij}(\rho_{i},\rho_{j})\Big{)}}_{i,j=1}^{n}\in{\mathbb{R}}^{n\times n}.$$

$${\mathbf{f}}_i$$

$$\sum _{j=1}^n{\tau }_{ij}$$

$${\mathbf{f}}_i=\sum _{j=1}^n{\tau }_{ij}({\mathbf{u}}_i-{\mathbf{u}}_j),i=1,...,n.$$

$$\zeta =-\sum _{i=1}^n{\mathbf{u}}_i\cdot {\mathbf{f}}_i=-\frac{1}{2}\sum _{i,j=1}^n{\tau }_{ij}|{\mathbf{u}}_i-{\mathbf{u}}_j{|}^2.$$

$$\zeta =\frac{1}{2}\sum _{i,j=1}^n{\lambda }_{ij}({\rho }_i,{\rho }_j){\rho }_i{\rho }_j|{\mathbf{v}}_i-{\mathbf{v}}_j{|}^2.$$

$${\mathbf{f}}_i=-\sum _{j=1}^n{\lambda }_{ij}({\rho }_i,{\rho }_j){\rho }_i{\rho }_j({\mathbf{v}}_i-{\mathbf{v}}_j).$$

$${\partial }_t({\rho }_i{\mathbf{v}}_i)+\mathrm{div}({\rho }_i{\mathbf{v}}_i{\mathbf{v}}_i^{\top }+p_i({\rho }_i){\mathcal{I}}_d)=-\sum _{j=1}^n{\lambda }_{ij}({\rho }_i,{\rho }_j){\rho }_i{\rho }_j({\mathbf{v}}_i-{\mathbf{v}}_j).$$

$$-\sum _{j=1}^n{\lambda }_{ij}({\rho }_i,{\rho }_j){\rho }_i{\rho }_j({\mathbf{v}}_i-{\mathbf{v}}_j)-{\lambda }_{i,\mathrm{pm}}({\rho }_i,{\rho }_{\mathrm{pm}}){\rho }_{\mathrm{pm}}{\rho }_i({\mathbf{v}}_i-{\mathbf{v}}_{\mathrm{pm}}).$$

$$M_i=M_i({\rho }_{\mathrm{pm}}):={\lambda }_{i,\mathrm{pm}}({\rho }_{\mathrm{pm}}){\rho }_{\mathrm{pm}}.$$

Full text

For inviscid flow the model turns out to be a nonlinear system of hyperbolic balance laws. We show that the dissipative structure of the Maxwell-Stefan operator permits to guarentee the existence of global classical solutions for initial data close to equilibria. Furthermore it is proven by relative entropy techniques that solutions of the Darcy-scale model tend in a certain long-time regime to solutions of a parabolic limit system.

AMS subject classifications: 35L65,76N10,76S05

1 Introduction

Multi-component flows in porous media appear in various fields of applications such as fuel cells, oxygen sensors, and respiratory airways [1, 3, 11, 26]. To highlight the modeling challenge let us focus on the last example. As stated in [3], the bronchial tree can be divided into two parts. In the lower part the velocity of the air is very small, such that the dynamics of the gas mixture is mainly dictated by diffusive effects. For the treatment of certain diseases of the lung, a gas mixture (Heliox) is used to improve the patient’s well-being. Mathematical models can be used to analyze how to achieve the greatest benefit for the patient. In this situation the classical Fickian diffusion law is too simplistic. Important effects, for instance uphill diffusion ([15]), cannot be covered by this approach. By uphill diffusion we mean flux from regions of low concentrations to ones with high concentration, see [4] and references therein. Let us additionally note that Duncan & Toor have given an experimental example of a three-component gas mixture in [9], which clearly demonstrates the uphill diffusion effect. A generalization of the Fickian approach is needed which roots in the classical works of Maxwell [20] and Stefan [22]. It has led to the concept which is nowadays called Maxwell–Stefan diffusion exploiting binary interactions between different species of the mixture. This approach captures more complex diffusive effects, but leads to a coupled nonlinear system of partial differential equations and is therefore mathematically more challenging.

In this paper we provide a mathematical model for compressible multi-component flow in fully saturated porous media on the Darcy-scale. This model takes the form of a nonlinear hyperbolic balance law, therefore classical solutions might fail to exist globally [7]. However, we show that dissipative effects of the Maxwell–Stefan diffusion and porous media friction suffice to ensure the classical wellposedness for initial data close to equilibrium. Note that we account for the effect of the solid skeleton in the porous medium like in dusty gas models from [19, 24]. It is regarded as an additional component of the mixture with zero velocity and constant density. In contrast to these dusty gas models, which rely on the kinetic theory of gases, we use the continuum thermodynamics framework as developed in [2]. The solutions of the resulting system automatically satisfy an entropy condition and hence the second law of thermodynamics. If we set the Maxwell–Stefan coupling terms to zero, we obtain a system of uncoupled equations which correspond to compressible Euler equations with friction. This model is used for single-component flow through porous media. It has been shown in [10, 12, 13, 16, 17] that the solutions to this system tend in a long-time limit to the solution of a parabolic porous media equation. We establish a corresponding result for the multi-component case that results in a parabolic system of porous media equations. This system is similar to the multi-component system of [14].

The paper is organized as follows. In Section 2 we derive the governing equations for compressible multi-component flow in porous media exploiting the continuum thermodynamics framework. We start with multiple Euler equations with friction, which are coupled by a right hand side using Maxwell–Stefan cross-diffusion terms. In the following Section 3 we investigate the existence of smooth solutions to this system. The dissipative effects due to Maxwell–Stefan diffusion and friction are shown to fit exactly into the wellposedness theory of [25] for general hyperbolic balance laws (see Theorem 1). The second major goal of the paper on the existence of a parabolic limit system is pursued in Section 4. We use a relative entropy framework to prove our main Theorem 5. For this purpose we adapt techniques from [16], where the convergence of the compressible Euler system with friction to the porous media equation has already been proven. We finally conclude in Section 5.

2 Compressible Flow in Porous Media

We first review in Section 2.1 a single-component model and highlight available analytical results which we will generalize to the multi-component case in the remainder of the paper. The multi-component model itself is derived using fundamental thermodynamical principles in Section 2.2. It follows the presentation in [2] for free flow problems extending the classical Fickian diffusion modeling to a Maxwell–Stefan approach.

2.1 Single-Component Flow

Compressible single-component flow with friction induced by the resistance of the solid skeleton in a porous medium can be described on a macroscopic averaged scale by the Euler–Darcy model, see [18] and references therein. We consider the isothermal situation with constant temperature. For $t\in(0,T)$, $T>0$, and $\mathbf{x}$ from the entire ${\mathbb{R}}^{d}$ the unknowns mass density $\rho=\rho(\mathbf{x},t)>0$ and momentum $\mathbf{m}=\rho\mathbf{v}\in{\mathbb{R}}^{d}$ with velocity $\mathbf{v}=\mathbf{v}(\mathbf{x},t)\in{\mathbb{R}}^{d}$ satisfy in this case the system

[TABLE]

Here $p=p(\rho)$ is the smooth pressure function, $M>0$ is the mobility constant, and $\mathcal{I}_{d}$ is the $d$-dimensional unit matrix. Note that we use in (2.1) the same notation for the divergence operator as applied to vector- or matrix-valued functions, see Appendix A.

It is well-known that (2.1) is a hyperbolic system of nonlinear balance laws as long as the pressure is monotonically increasing. Shock-type singularities might evolve in finite time regardless of the initial data’s regularity. It is of overall importance that (2.1) is endowed with a entropy-entropy flux pair (see (2.2) below) which can be used to ensure an appropriate form of the second law of thermodynamics for classical as well as weak solutions. Precisely, solutions of (2.1) are required to satisfy the entropy inequality

[TABLE]

for the entropy-entropy flux pair $(\eta,\mathbf{q})$ given by

[TABLE]

Here the free energy density $\rho\psi$ is determined from the Gibbs–Duhem relation

[TABLE]

While the nonlinear flux in (2.1) can drive shock waves, the dissipative effect of the friction term might suffice to counteract the destabilizing effect of the flux. Depending on the initial data and the size of $M$ the initial value problem (IVP) for (2.1) can have in fact global smooth solutions (see e.g. [10]). We will show that a similar result holds for the multi-component case. Furthermore, the the dissipative friction effect leads to certain limit regimes such that (2.1) changes type in the limit. We consider a long-time and large-mobility regime in (2.1), i.e., the time $t$ is scaled by a small parameter $\varepsilon>0$ and the mobility $M$ by $\varepsilon^{-1}$. After rescaling (2.1) and renaming the variables in an obvious way (2.1) is recasted in the form

[TABLE]

In [16] (but see also [13]) it has been shown that the sequence of densities ${\rho}_{\varepsilon}>0$ solving the initial value problems for (2.4) converge for $\varepsilon\to 0$ towards a solution $\bar{\rho}$ of the IVP for the porous medium equation

[TABLE]

In other words the hyperbolic balance laws turn into a parabolic evolution with much more regular solution behavior. We aim at a corresponding result for the multi-component case (for the parabolic system (2.26) that reduces to (2.5) in the one-component case).

2.2 Multi-Component Flow

While single-component flow in a porous medium is well understood, much less is known for multi-component flow. As long as bulk viscosity is neglected standard model approaches take the form of the Euler equations with a damping term in the momentum equations like (2.1). However, in the case of porous media and multi-component gaseous mixtures, inter-component viscosity effects become important which do not occur in the single-component case. The classical Fickian approach does not suffice to describe these diffusion phenomena. As a possible remedy we favor in this paper a Maxwell–Stefan ansatz and, in order to derive governing equations in a thermodynamically consistent way, follow the work of Bothe & Dreyer [2] for free multi-component flow.

2.2.1 Multi-Component Flow and Maxwell–Stefan Diffusion

Let a fluid mixture consist of $n\in{\mathbb{N}}$ components $A_{1},...,A_{n}$ with corresponding mass densities $\rho_{i}=\rho_{i}(\mathbf{x},t)>0$ and velocities $\mathbf{v}_{i}=\mathbf{v}_{i}(\mathbf{x},t)\in{\mathbb{R}}^{d},$ $i=1,...,n$.

We define the total mass density $\rho$ and the barycentric velocity $\mathbf{v}$ (not to be interchanged with the single-component case in Section 2.1) as

[TABLE]

Further, we define the diffusion velocities

[TABLE]

We ignore mass exchange as well as exterior forces. Restricting ourselves to the case of a simple mixture, the component pressures $p_{i}$ depend on $\rho_{i}$ only, i.e. they satisfy $p_{i}=p_{i}(\rho_{i})$. For $i=1,\dots,n$ we start then from the partial balances of mass and momentum given by

[TABLE]

Here $\mathbf{f}_{i}\in{\mathbb{R}}^{d}$ states the momentum production due to diffusive mixing, later to be specified with the Maxwell–Stefan ansatz. As a natural requirement the conservation law for total momentum has to hold, which implies the condition

[TABLE]

The crucial part is now to find an expression for $\mathbf{f}_{i}$ such that with the (physical) entropy production $\zeta$ (see (2.8) below) the second law of thermodynamics holds true.

We introduce for each component $A_{i}$ a strictly convex free energy density $h_{i}(\rho_{i})=\rho_{i}\psi_{i}(\rho_{i})$ that relates to the partial pressure $p_{i}(\rho_{i})$ via the Gibbs–Duhem equations (see (2.3) for the single-component velocity)

[TABLE]

Thus, the strict convexity of $\rho_{i}\psi_{i}$ implies

[TABLE]

Moreover, the function

[TABLE]

acts as mixture free energy for simple mixtures.

For the special case of simple isothermal inviscid fluid mixtures without chemical reactions the entropy production $\zeta$ of some solution $(\rho_{1},\ldots,\rho_{n},\mathbf{m}_{1}^{\top},\ldots,\mathbf{m}_{n}^{\top})$ of (2.6a), (2.6b) is derived in [2] and reads as

[TABLE]

The nonnegativity of the entropy production $\zeta$ leads then to the requirement

[TABLE]

In the following we make the Maxwell–Stefan ansatz for $\mathbf{f}_{i}$ to guarantee that (2.9) holds true. Let

[TABLE]

be a positive-definite matrix. With (2.10) we set

[TABLE]

In order to make the right hand side in (2.11) symmetric regarding the components, we extend $\tilde{\mathcal{T}}\in{\mathbb{R}}^{(n-1)\times(n-1)}$ to the Maxwell–Stefan matrix $\mathcal{T}:=(\tau_{ij})_{i,j=1}^{n}\in{\mathbb{R}}^{n\times n}$ [23] by

[TABLE]

Additionally let

[TABLE]

In the case of binary interactions the matrix $\mathcal{T}$ is symmetric and can be shown to be positive semi-definite provided (2.10) holds, see [2].

The following ansatz for the components of $\tilde{\mathcal{T}}$ is made to match the requirements (2.12), (2.13). For some numbers $\lambda_{ij}\geq 0$ with $\lambda_{ij}(\rho_{i},\rho_{j})=\lambda_{ji}(\rho_{j},\rho_{i})$ for $i\neq j$ we define

[TABLE]

With $\lambda_{ii}=-\sum_{j=1,j\neq i}^{n}\lambda_{ij}\frac{\rho_{j}}{\rho_{i}}$ we introduce the negative semi-definite matrix

[TABLE]

With the definitions (2.12) we infer from (2.11) the relations

[TABLE]

Thus, by replacing $\mathbf{u}_{n}$ with $\mathbf{u}_{i}$ in (2.16), we obtain a symmetrical version of (2.11), namely

[TABLE]

With (2.17) and the symmetry of $\mathcal{T}$ the entropy production (2.8) can be written as

[TABLE]

Obviously, condition (2.13) is necessary to achieve $\zeta\geq 0$. Due to (2.14) the entropy production reads as

[TABLE]

Finally, with the Maxwell–Stefan ansatz the constitutive law for the momentum production $\mathbf{f}_{i}$ results from (2.14) and (2.17) as

[TABLE]

Note that in (2.19) the diffusion velocities are replaced by the velocities of the corresponding component.

With this result the partial momentum balances (2.6b) read as

[TABLE]

2.2.2 Compressible Multi-Component Flow in Porous Media

So far we considered a free flow problem. We realize the porous medium matrix as a static component $A_{\mathrm{pm}}$ of the mixture with velocity $\mathbf{v}_{\mathrm{pm}}=0$ and density $\rho_{\mathrm{pm}}=\mathrm{const.}$

The component $A_{\mathrm{pm}}$ needs no equations for the mass and momentum balance. However, we have to take into account the effects on the other components. Hence, the sum from (2.19) extends to

[TABLE]

In the sequel we ignore the explicit dependence of $\lambda_{i,\text{pm}}$ on the component densities $\rho_{i}$ and proceed with the mobility constants

[TABLE]

Then we arrive at our final system which reduces in the single-component case to (2.1).

Define the density vector $\mathbf{r}=(\rho_{1},\ldots,\rho_{n})^{\top}$ and with $\mathbf{m}_{i}=\rho_{i}\mathbf{v}_{i}$ for $i=1,\ldots,n$ the momentum vector $\mathbf{m}=(\mathbf{m}_{1}^{\top},\ldots,\mathbf{m}_{n}^{\top})^{\top}$. We search for the function $U=(\mathbf{r}^{\top},\mathbf{m}^{\top})^{\top}$ with values in the state space

[TABLE]

that satisfies the multi-component Euler–Darcy system with Maxwell–Stefan type diffusion

[TABLE]

subject to the initial condition

[TABLE]

Due to the arguments from Section 2.2.1 we observe that the second law of thermodynamics is automatically satisfied along smooth solution trajectories of (2.22), (2.23).

Remark \thethm.

When including the porous medium, the condition (2.16b) contains the summand of the porous medium part as well. We neglect this fact. That means precisely that we neglect the effect of the porous medium on the conservation of momentum (2.19). This is in accordance with the single-component case. The effect of the porous medium has the character of a body force. 2. 2.

It turns out in experiments that the $\lambda_{ij}$ are only weakly dependent on the mixture. Often affine functions suffice to describe this dependence, see [23]. 3. 3.

The structure of the porous medium is only captured in the scalar parameter $M_{i}$. If the porous medium is not homogeneous and isotropic, one should allow for matrix-valued parameters with spatial dependence. 4. 4.

The terms $M_{i}$ scale with the density of the porous medium, which is significantly larger than the densities of a gaseous mixture. Hence, typically it holds $M_{i}\gg\lambda_{ij}$.

Similar like in Section 2.1 we consider a long-time/large-mobility/large-diffusion regime for (2.22). To be precise, let $\bar{x},\,\bar{t},\,\bar{\rho}>0,\,\bar{v},\,\bar{p},\,\bar{M},\,$ and $\bar{\lambda}$ be the characteristic scales of the corresponding quantities. The long-time/large-mobility/large-diffusion regime is now obtained from

[TABLE]

using a small parameter $\varepsilon>0$. After rescaling the system $\eqref{Euler-MS}$ and renaming the unknown as $U^{\varepsilon}=(\mathbf{r}^{\varepsilon\top}=(\rho_{1}^{\varepsilon},\ldots,\rho_{n}^{\varepsilon}),\mathbf{m}^{\varepsilon\top}=(\mathbf{m}_{1}^{\varepsilon\top},\ldots,\mathbf{m}_{n}^{\varepsilon\top}))^{\top}$ it reads in this regime as

[TABLE]

with the ($\varepsilon$-dependent) initial conditions

[TABLE]

As for the single-component case we will show that the multi-component case admits global smooth solutions exploiting the dissipative effect due to friction and the Maxwell–Stefan diffusion. The other major goal is to prove that the density component sequence ${\{\mathbf{r}^{\varepsilon}\}}_{\varepsilon>0}$ of solutions of the IVP (2.24), (2.25) converges for $\varepsilon\to 0$ to the vector-valued density field $\mathbf{\bar{r}}=(\bar{\rho}_{1},...,\bar{\rho}_{n})^{\top}$ solving the system of porous medium equations

[TABLE]

subject to the initial conditions

[TABLE]

In (2.26) we used the vector-valued pressure

[TABLE]

and the matrices

[TABLE]

The matrix $\Lambda$ from (2.15) is negative semi-definite. This implies that $\mathcal{\tilde{B}}$ is positive definite due to $M_{i}>0$ and the positivity of the densities. In particular $\mathcal{B}$ turns then out to be positive definite as the Kronecker product of two positive definite matrices. For the definition of the generalized gradient/divergence operators in the system (2.26) and the Kronecker matrix product $\otimes$ in (2.29) we refer to Appendix A.

Remark \thethm.

For the single component case $n=1$ the system (2.26) reduces to the porous media equation (2.5). 2. 2.

If no porous medium is present, that is $M_{i}=0$, the system (2.26) in this framework corresponds for perfect gas laws to the following version of the Maxwell–Stefan equations formulated for the molar concentrations $c_{i}$ often seen in the literature, e.g. in[14]:

[TABLE]

Here $D_{ij}=\dfrac{R}{c\mathcal{M}_{i}\mathcal{M}_{j}\lambda_{ij}}$, with the ideal gas constant $R$, total molar concentration $c=\sum_{i=1}^{n}c_{i}$ and molar masses $\mathcal{M}_{i}$.

3 Existence of Smooth Solutions in Multiple Space Dimensions

The main result in this section is Theorem 1 on the classical wellposedness of the IVP for system (2.22). To this end we propose an entropy concept for (2.22) and adapt a result of Yong [25] on hyperbolic balance laws, see Appendix B. It exploits dissipative effects of the balance terms that counteract the development of singularities driven by the hyperbolic flux [7]. To state and prove our main result Theorem 1 below we summarize all assumptions on the system (2.22) according to the notations from Section 2.

Assumption \thethm.

• (i)

The functions $\lambda_{ij}\in C^{\infty}((0,\infty)^{2},{\mathbb{R}}),\ i,j=1,\dots,n$, satisfy (2.14) and (2.16b).

• (ii)

The symmetric matrix \Lambda(\mathbf{r})={\Big{(}\lambda_{ij}(\rho_{i},\rho_{j})\Big{)}}_{i,j=1}^{n} is negative semi-definite for all $\mathbf{r}\in(0,\infty)^{n}$.

• (iii)

The free energy densities $h_{i}=\rho_{i}\psi_{i}\in C^{3}((0,\infty))$ are strictly convex for $i=1,\dots,n$.

• (iv)

The mobility constants $M_{i}$ are positive for $i=1,\ldots,n$.

Theorem 8 applies to general hyperbolic balance laws. With $U=(\rho_{1},\dots,\rho_{n},(\mathbf{m}_{1}^{\top},\dots,\mathbf{m}_{n}^{\top}))^{\top}$ we can rewrite (2.22) in this form, that is

[TABLE]

The fluxes $F_{\alpha}(U)\in{\mathbb{R}}^{(d+1)n}$ and the source $\mathbf{s}(\mathbf{r},\mathbf{m})\in{\mathbb{R}}^{nd}$ are given by

[TABLE]

Here we used $\mathbf{m}_{i}=(\mathbf{m_{i}}^{(1)},\ldots,\mathbf{m_{i}}^{(d)})^{\top}$ and $\mathbf{e}_{\alpha}$ denotes the $\alpha$-th unit vector. Furthermore, an entropy-entropy flux pair $(\eta,\mathbf{q})\in C^{2}(G)$ for (3.1) on the state space $G$ from (2.21) is required. Following [7] the tuple $(\eta,\mathbf{q})$ is called an entropy-entropy flux pair to the system (2.22) provided $\operatorname{D}^{2}\eta(U)$ is positive-definite and the compatibility conditions

[TABLE]

are satisfied for all $U\in G$. Motivated by the considerations in Section 2 we suggest for (2.22) the functions

[TABLE]

Note that $\eta$ in (3.4) is obviously strictly convex, by Assumption 3 (iii).

Theorem 1** (Global classical wellposedness of the IVP for (2.22)).**

Let $s\geq\lfloor d/2\rfloor+2$ and let Assumption 3 hold. Consider a static equilibrium solution $\hat{U}\in G$ of (2.22) of the form

[TABLE]

Then there exists a constant $c_{1}>0,$ such that for all $U_{0}\in H^{s}({\mathbb{R}}^{d})$ with

[TABLE]

and all $T>0$ the IVP (2.22), (2.23) has a unique solution $U\in C([0,T),H^{s}({\mathbb{R}}^{d}))$ taking values in the state space $G$.

Additionally the solution $U$ satisfies the entropy inequality

[TABLE]

with $\zeta=\sum_{i,j=1}^{n}\frac{\lambda_{ij}(\rho_{i},\rho_{j})}{2\rho_{i}\rho_{j}}|\rho_{j}\mathbf{m}_{i}-\rho_{i}\mathbf{m}_{j}|^{2}$ (see (2.18)).

Proof 3.1**.**

For better readability we omit the argument $(\rho_{i},\rho_{j})$ of the functions $\lambda_{ij}$. With (3.1) we meet the setting (B.2) to apply Theorem 8. We compute for the equilibrium state $\hat{U}$ from (3.5) the Jacobian

[TABLE]

From Assumption 3 (i) we have $\hat{\rho}_{i}\lambda_{ii}=-\sum_{j=1,j\neq i}^{n}\hat{\rho}_{j}\lambda_{ij}$ which implies

[TABLE]

with $\mathcal{M}=\operatorname{diag}(M_{i})\in{\mathbb{R}}^{n\times n}$, $\mathcal{R}=\operatorname{diag}(\hat{\rho}_{i})\in{\mathbb{R}}^{n\times n}$ being positive-definite. The Jacobian $\operatorname{D}_{\mathbf{m}}\mathbf{s}(\hat{U})\prec 0$ is in particular a regular matrix which implies condition 1. in Theorem 8.

For the entropy $\eta$ from (3.4) and any open set $\hat{\mathcal{G}}\subset G$ containing $\hat{U}$ we have $\operatorname{D}\eta(U)\in{\mathbb{R}}^{n(d+1)},\operatorname{D}^{2}\eta\in{\mathbb{R}}^{n(d+1)\times n(d+1)}$, $U\in\hat{\mathcal{G}}$. The derivatives of $\eta$ read with (2.7) as

[TABLE]

Fro the definition of the operator $\operatorname{blockdiag}$ and the Kronecker product $\otimes$ we refer to the Appendix B. From the definition of the fluxes $F_{\alpha}$ in (3.2) we compute

[TABLE]

where $\mathbf{e}_{\alpha}$ denotes the $\alpha$-th unit vector and $1_{d}:=(1,\dots,1)^{\top}\in{\mathbb{R}}^{d}$. We see that the matrix

[TABLE]

*is symmetric. Altogether we have verified condition 2. in Theorem 8, up to now for any open set $\hat{\mathcal{G}}$ that contains the equilibrium state $\hat{U}$.

To check the third condition 8 let $U,\hat{U}\in\hat{\mathcal{G}}$. We use the symmetry of $\Lambda$ and (2.18) to obtain*

[TABLE]

Furthermore we now choose $\mathcal{G}$ as a compact, convex subset of $\hat{\mathcal{G}}$ such that we have for all $U\in\mathcal{G}$ the estimate

[TABLE]

Hence, we get from (3.7) with $c_{\mathcal{G}}^{-1}=2\max\left\{{\hat{c}_{\mathcal{G}}},{\tilde{c}_{\mathcal{G}}}\right\}$ the inequality

[TABLE]

*which implies the third condition.

Finally,*

[TABLE]

with $\mathcal{A}(\hat{U})=0_{nd\times n}$. We obtain

[TABLE]

The lower right block of this matrix is invertible as shown above. Consequently,

[TABLE]

Due to the zero block in $\operatorname{D}F_{j}(\hat{U})$, the corresponding eigenvectors must have non-zero entries at the $n+1$-th to $n(d+1)$-th position. Therefore, the last condition of Theorem 8 holds.

We verified all the conditions of Theorem 8. Hence, the system (2.22) with $U_{0}$ as initial value has a unique solution $U=U(\mathbf{x},t)\in C([0,T),H^{s}({\mathbb{R}}^{d})),\,s\geq\lfloor d/2\rfloor+2$. The entropy inequality (3.6) is a consequence of (3.3), (3.7), and $\operatorname{D}\eta(\hat{U})S(U)=0$.

Remark 2**.**

Note that due to $s\geq\lfloor d/2\rfloor+2$ we have with the Sobolev embedding theorem even $U\in C^{1}({\mathbb{R}}^{d}\times(0,T))$. 2. 2.

The symmetry of $\operatorname{D}^{2}\eta(U)\operatorname{D}F_{j}(U)$ follows directly from the compatibility condition (3.3) of the entropy-entropy flux pair and the strict convexity of $\eta$. However, since the matrices are needed in the proof anyway, we checked this property by hand.

4 Convergence to the Parabolic Limit System

The goal of this section is to prove the convergence of solutions of (2.24), (2.25) to solutions of an IVP for the parabolic limit system (2.26) as $\varepsilon$ tends to zero. Due to Theorem 1 there exists for each $\varepsilon>0$ a unique global solution ${U}^{\varepsilon}=(\mathbf{r}^{\varepsilon\top},\mathbf{m}^{\varepsilon\top})^{T}$ to the IVP for the $\varepsilon$-scaled system (2.24). However, the convex set $\mathcal{G}$ might depend on $\varepsilon$ such that the set of admissible initial conditions could shrink to the equilibrium for $\varepsilon\to 0$. The techniques of [21] allow to show that there exists a time interval independent of $\varepsilon$ with solutions $U^{\varepsilon}$ existing. Hence, we assume that there is a time $T>0$ and a compact, convex set $\mathcal{G}$ such that for all $\varepsilon>0$ the solutions $U^{\varepsilon}$ exist on the interval $(0,T)$ and are contained in $\mathcal{G}$.

Our convergence proof relies on the relative entropy method which goes back to [5, 6] and [8]. This technique only requires one solution to be a strong (in fact Lipschitz continuous) solution, whereas the other can be a discontinuous entropy solution. We regard the solutions to (2.24) as weak solutions and the solution of the limit system as strong solution. Here we rely on the technical framework that has been established in [16]. We start to prove a dissipation relation (Proposition 3) for so-called relative entropies in Section 4.1 and conclude the convergence estimate with the main result in Theorem 5 of Section 4.2.

In the following we omit again the arguments in $\lambda_{ij}$. With a slight misuse of notation we use the expression $(\mathbf{r}^{\varepsilon},\mathbf{m}^{\varepsilon})$ for the solution $U^{\varepsilon}$.

4.1 The Relative Entropy Estimate

Let us consider (2.24), (2.25) for $\varepsilon>0$. To obtain a convergence estimate for the solutions $(\mathbf{r}^{\varepsilon},\mathbf{m}^{\varepsilon})$ of (2.24), (2.25) we start to fix well prepared functions for the initial conditions in (2.25) and (2.27) on the entire ${\mathbb{R}}^{d}$. Let a number $R_{0}>0$ and $\hat{\mathbf{r}}\in(0,\infty)^{n}$ be given. We restrict the initial datum $(\mathbf{r}^{\varepsilon}_{0},\mathbf{m}^{\varepsilon}_{0})\in L^{\infty}({\mathbb{R}}^{d})$ in (2.25) to take values in $\mathcal{G}$ and to satisfy

[TABLE]

For the initial datum $\mathbf{\bar{r}}_{0}\in C^{3}({\mathbb{R}}^{d})$ of the limit equation (2.26) we impose the analogous condition

[TABLE]

Using the entropy-entropy flux pair $(\eta,\mathbf{q})$ from (3.4) we define an entropy solution $(\mathbf{r}^{\varepsilon},\mathbf{m}^{\varepsilon})\in L^{\infty}({\mathbb{R}}^{d}\times(0,T))$ of (2.24), (2.25) as a weak solution of (2.24), (2.25) that takes values in $\mathcal{G}$ and satisfies

[TABLE]

in $D^{\prime}({\mathbb{R}}^{d}\times[0,T))$. Note that entropy flux scales with $\varepsilon^{-1}$ according to the flux scaling in (2.24).

Further, let $\mathbf{\bar{r}}\in C^{3,1}({\mathbb{R}}^{d}\times(0,T))\coloneqq\{\mathbf{g}\ |\ \mathbf{g}(\cdot,t)\in C^{3}({\mathbb{R}}^{d}),t\in(0,T),\ \mathbf{g}(\mathbf{x},\cdot)\in C^{1}((0,T)),\mathbf{x}\in{\mathbb{R}}^{d}\}$ with $\mathbf{\bar{r}}\in(0,\infty)^{n}$ be a smooth solution of (2.26) (2.27). We observe that $\mathbf{\bar{r}}$ satisfies for all $\varepsilon>0$ the expanded but equivalent system

[TABLE]

Recall that the matrix $\mathcal{B}$ has been defined in (2.29). The regularity of $\mathbf{\bar{r}}$ implies $\mathbf{\bar{m}}\in C^{2,1}({\mathbb{R}}^{d}\times(0,T))$ for the momentum. Note that $\mathbf{\bar{m}}$ depends on $\varepsilon$ which is suppressed in the notation. The equivalent formulation (4.4) of system (2.26) induces for the evaluation of $\mathbf{\bar{m}}$ at zero time by (4.2) the compatibility condition

[TABLE]

For this choice of $\mathbf{\bar{m}}$ we define now the relative entropy expression

[TABLE]

with

[TABLE]

The relative entropy flux is defined by

[TABLE]

with $F$ being a vectorial collection of the momentum fluxes given by

[TABLE]

In the last formula we made use of the notation

[TABLE]

which appears frequently in the sequel, where also Lemma 7 will be used often.

After artificially expanding the system (2.26) to obtain (4.12), we are able to compare the solutions $(\mathbf{r}^{\varepsilon},\mathbf{m}^{\varepsilon})$ and $(\bar{\mathbf{r}},\bar{\mathbf{m}})$ of (2.24) and (4.12), respectively.

Proposition 3**.**

Let Assumption 3 hold, let the pressure $p_{i}$ satisfy $\eqref{passum1}$, and let the initial functions $(\mathbf{r}^{\varepsilon}_{0},\mathbf{m}^{\varepsilon}_{0})\in L^{\infty}({\mathbb{R}}^{d})$ and $\mathbf{\bar{r}}_{0}\in C^{3}({\mathbb{R}}^{d})$ satisfy (4.1), (4.2).

*Consider for $\varepsilon>0$ an entropy solution $(\mathbf{r}^{\varepsilon},\mathbf{m}^{\varepsilon})\in L^{\infty}({\mathbb{R}}^{d}\times(0,T))$ of (2.24), (2.25) and a smooth solution $(\bar{\mathbf{r}},\bar{\mathbf{m}})\in C^{3,1}({\mathbb{R}}^{d}\times[0,T))\times C^{2,1}{\mathbb{R}}^{d}\times[0,T))$ of (2.26), (2.27), supposed to take values in a convex, compact set $\mathcal{G}\subset G$.

Then we have the estimate

[TABLE]

with

[TABLE]

Before we present the proof of Proposition 3 we summarize two remarks on the scaling of the remainder terms $Q_{\varepsilon}$ and $E_{\varepsilon}$ with respect to $\varepsilon$ which will be needed in Section 4.2.

Remark 4**.**

The first factor of $Q_{\varepsilon}$ depends only on $(\bar{\mathbf{r}},\bar{\mathbf{m}})$. Although $\bar{\mathbf{m}}$ involves $\varepsilon$ the factor is independent of $\varepsilon$, i. e.,

[TABLE] 2. 2.

Recalling the definition of $\tilde{B}$ from (2.29) the smoothness of $\bar{\mathbf{r}}$ and $\mathbf{p}$ implies

[TABLE]

Hence the vector $\bar{\mathbf{e}}_{\varepsilon}=O(\varepsilon)$ in $E_{\varepsilon}$ is of order $O(\varepsilon)$.

Proof (of Proposition 3). To simplify notations we may omit the index $\varepsilon$ and write $(\mathbf{r},\mathbf{m})=(\mathbf{r}^{\varepsilon},\mathbf{m}^{\varepsilon})$. The entropy solution $(\mathbf{r},\mathbf{m})$ of the IVP for (2.24) satisfies the inequality (4.3).

In order to derive a similar expression for the solution $\mathbf{\bar{r}}$ of the IVP for (2.26) we rewrite the equivalent system (4.4) for the pairing $(\mathbf{\bar{r}},\mathbf{\bar{m}})$ further.

With $\lambda_{ii}r_{i}=-\sum_{i\neq j}\lambda_{ij}r_{j}$, (2.14), and (2.16b), it is easy to see that the solution $(\mathbf{\bar{r}},\mathbf{\bar{m}})$ of (4.4) also satisfies

[TABLE]

with $\bar{\mathbf{e}}_{\varepsilon}$ from (4.10) and $F$ from (4.8).

With (4.2) we see that $(\mathbf{\bar{r}},\mathbf{\bar{m}})$ satisfies in the sense of distributions

[TABLE]

Before we use the entropy relations (4.3) and (4.13) we return to the weak formulations: We subtract the weak formulations of (4.12) from the weak formulation for (2.24) and obtain for the mass balance equations

[TABLE]

Using the definition of $\mathbf{\bar{e}}_{\varepsilon}$ in (4.10) yields for the momentum components

[TABLE]

Here, $\boldsymbol{\varphi}$ and $\boldsymbol{\theta}$ are vector-valued test functions with compact support in ${\mathbb{R}}^{d}\times[0,T)$. We make with some function $\psi\in C^{\infty}_{0}({\mathbb{R}}^{d}\times[0,T))$ the specific choices

[TABLE]

which lead in (4.14) and (4.15) to

[TABLE]

The term $J_{\varepsilon}$ in (4.16) is defined as

[TABLE]

Combining the entropy inequality (4.3) for $(\mathbf{r}^{\varepsilon},\mathbf{m}^{\varepsilon})$, the entropy equation (4.13) for $(\bar{\mathbf{r}},\mathbf{\mathbf{m}})$ and the relation (4.16) the definition of the relative entropy-entropy flux pair in (4.6), (4.1) implies that the inequality

[TABLE]

holds in the weak sense.

In the term $J_{\varepsilon}$ we use (4.13) and Lemma 7 to compute the time derivative of $\nabla_{\mathbf{m}}\eta(\mathbf{\bar{r}},\mathbf{\bar{m}})$ by the chain rule which leads to

[TABLE]

Finally we proceed with this expression for $J_{\epsilon}$ and deduce

[TABLE]

The last line follows from the definitions in (4.10) and concludes the proof.

4.2 The Convergence Estimate

In this section we make an additional assumption on the pressure. Let there exist constants $a_{i}>0,\ i=1,\dots,n$, such that

[TABLE]

The condition (4.17) is satisfied for e.g. the isentropic pressure laws $p_{i}(r)=k_{i}r^{\gamma_{i}}$ ($\gamma_{i}\geq 1,\ k_{i}>0$) with any choice of $a_{i}>0$.

Due to (4.17) we have

[TABLE]

Note that $\displaystyle p_{i}(\rho_{i}|\bar{\rho}_{i})=p_{i}(\rho_{i})-p_{i}(\bar{\rho}_{i})-p_{i}^{\prime}(\bar{\rho}_{i})(\rho_{i}-\bar{\rho}_{i})=(\rho_{i}-\bar{\rho}_{i})^{2}\int_{0}^{1}\int_{0}^{\tau}p^{\prime\prime}(s\rho_{i}+(1-s)\bar{\rho}_{i})\operatorname{~{}d\!}s\operatorname{~{}d\!}\tau.$

Hence with

[TABLE]

the inequality (4.18) implies with some $c>0$

[TABLE]

with $F(\mathbf{r},\mathbf{m}|\mathbf{\bar{r}},\mathbf{\bar{m}})=F(\mathbf{r},\mathbf{m})-F(\mathbf{\bar{r}},\mathbf{\bar{m}})-\operatorname{D}F(\mathbf{\bar{r}},\mathbf{\bar{m}})(\mathbf{r}-\mathbf{\bar{r}},\mathbf{m}-\mathbf{\bar{m}}).$

We need to introduce a slightly different entropy-entropy flux pair to get a convergence estimate that corresponds to convergence in standard Lebesgue spaces. Subtracting the constant $\eta(\mathbf{\hat{r}},\mathbf{0})=\sum_{i=1}^{n}h_{i}(\hat{\rho_{i}})$ from the entropy $\eta$ we obtain a modified entropy-entropy flux pair $(\tilde{\eta},\tilde{\mathbf{q}})$ with the property $\tilde{\eta}(\hat{\mathbf{r}},0)=0$ by setting

[TABLE]

Since (2.24) is a hyperbolic balance law, due to (4.1) and the uniform bound in $\mathcal{G}$ the functions $(\mathbf{r}^{\varepsilon}-\hat{\mathbf{r}},\mathbf{m}^{\varepsilon})$ have uniform compact support. Then again the uniform boundedness implies that there are constants $K_{1},K_{2}>0$ such that for any $\varepsilon>0$ the entropy solution $(\mathbf{r}^{\varepsilon},\mathbf{m}^{\varepsilon})$ of (2.24), (2.25) satisfies

[TABLE]

As discussed in the introduction to Section 4 we will consider a classical solution $\bar{\mathbf{r}}$ of (2.26). Let $(\bar{\mathbf{r}},\bar{\mathbf{m}})\in C^{3,1}({\mathbb{R}}^{d}\times[0,T))\times C^{2,1}{\mathbb{R}}^{d}\times[0,T))$ be a classical solution of (2.26) (respectively the equivalent system (4.4)), (2.27) with initial data satisfying (4.2). Since (2.26) is a regular parabolic system we can assume under corresponding conditions on $\bar{\mathbf{r}}_{0}$,that $(\bar{\mathbf{r}},\bar{\mathbf{m}})$ is also contained in $\mathcal{G}$. With the relative entropy $\tilde{\eta}(\mathbf{r}^{\varepsilon},\mathbf{m}^{\varepsilon}|\bar{\mathbf{r}},\bar{\mathbf{m}})={\eta}(\mathbf{r}^{\varepsilon},\mathbf{m}^{\varepsilon}|\bar{\mathbf{r}},\bar{\mathbf{m}})-\sum_{i=1}^{n}h_{i}(\hat{\rho_{i}})$, we measure the distance between the solutions $(\mathbf{r}^{\varepsilon},\mathbf{m}^{\varepsilon})$ and $(\bar{\mathbf{r}},\bar{\mathbf{m}})$ via the expression

[TABLE]

Note that the conditions (4.1), (4.2), and (4.5) show that $\varphi_{\varepsilon}(0)$ is finite. Due to the strict convexity of $\tilde{\eta}$ there are some constants $c,C>0$ (which depend on $\mathcal{G}$) such that

[TABLE]

holds for all vectors $(\mathbf{s},\mathbf{n}),(\bar{\mathbf{s}},\bar{\mathbf{n}})\in\mathcal{G}$. As a consequence (4.21) is compatible with the $L^{2}$-difference of the solutions $(\mathbf{r}^{\varepsilon},\mathbf{m}^{\varepsilon})$ and $(\bar{\mathbf{r}},\bar{\mathbf{m}})$. We can now state the final theorem.

Theorem 5** (Asymptotic behavior for (2.24)).**

Let Assumption 3 hold, let the pressure $p_{i}$ satisfy $\eqref{passum1}$, and let the initial functions $(\mathbf{r}^{\varepsilon}_{0},\mathbf{m}^{\varepsilon}_{0})\in L^{\infty}({\mathbb{R}}^{d})$ and $\mathbf{\bar{r}}_{0}\in C^{3}({\mathbb{R}}^{d})$ satisfy (4.1), (4.2).

*Consider for $\varepsilon>0$ an entropy solution $(\mathbf{r}^{\varepsilon},\mathbf{m}^{\varepsilon})\in L^{\infty}({\mathbb{R}}^{d}\times(0,T))$ of (2.24), (2.25) and a smooth solution $(\bar{\mathbf{r}},\bar{\mathbf{m}})\in C^{3,1}({\mathbb{R}}^{d}\times[0,T))\times C^{2,1}{\mathbb{R}}^{d}\times[0,T))$ of (2.26), (2.27), supposed to take values in a convex, compact set $\mathcal{G}\subset G$.

Then there exist constants $c_{i}>0,i=1,\dots,n$, such that for*

[TABLE]

we have the estimate

[TABLE]

Here $K>0$ is a constant that depends only on $T$, $\mathcal{G}$ and $\bar{\mathbf{r}}$ but not on $\varepsilon$.

Remark 6**.**

If the initial datum $(\mathbf{r}^{\varepsilon}_{0},\mathbf{m}^{\varepsilon}_{0})$ converges for $\varepsilon\to 0$ to $(\bar{\mathbf{r}}_{0},\bar{\mathbf{m}}(\cdot,0))$ in $L^{2}_{loc}({\mathbb{R}}^{d})$ the estimate (4.24) implies

[TABLE]

due to the compatibility relation (4.22). 2. 2.

The condition (4.23) holds especially for $\lambda_{ij}=0,\ i,j=1,\dots,n$ what corresponds exactly to **[16]**. In gaseous mixtures (4.23) is expected to hold, see Remark 2.2.2.

Proof (of Theorem 5). For the proof we write again $(\mathbf{r},\mathbf{m})=(\mathbf{r}^{\varepsilon},\mathbf{m}^{\varepsilon})$.

We consider the relative entropy statement from (4.9) in Proposition 3 which holds also for the entropy pair $(\tilde{\eta},\tilde{\mathbf{q}})$. As test function $\psi$ we make the choice $\psi(\mathbf{x},t)=\theta_{\kappa}(\tau)\omega_{R}(\mathbf{x})$ with $\theta_{\kappa}$ given for $\kappa>0$ by

[TABLE]

and with $\omega_{R}$ given for $R,\delta>0$ by

[TABLE]

By taking the limit $R\to\infty$, using the asymptotic properties (4.1), (4.2) of $(\mathbf{r},\mathbf{m})$ and $(\bar{\mathbf{r}},\bar{\mathbf{m}})$, and finally sending $\kappa\to 0$ we obtain using the definition of $\varphi_{\varepsilon}$ from (4.21) the inequality

[TABLE]

with $Q_{\varepsilon},E_{\varepsilon},$ and $R_{\varepsilon}$ from (4.10) in Proposition 3.

Due to Remark 4 and (4.19) it holds

[TABLE]

where $C_{1}>0$ depends on the $L^{\infty}$-norm of $\nabla\bar{\mathbf{r}}$ but not on $\varepsilon$. The error term $E_{\varepsilon}$ can be estimated for any number $C_{2}>0$ with Young’s inequality by

[TABLE]

Additionally, we have from (4.11) with $\bar{\mathbf{e}}_{\varepsilon}=O(\varepsilon)$ (see Remark 4) the inequality

[TABLE]

where the constant $C_{3}>0$ depends on $T,K_{1}$ from (4.20), $\mathcal{G}$, and also on $\bar{\mathbf{r}}$ through (4.11).

Plugging these estimates into (4.26) leads to

[TABLE]

We need the integral on the left hand side of the last estimate to be positive. The integrand reads as

[TABLE]

The term $R_{1,\varepsilon}$ is positive and scales with the mobilities $M_{i}$, whereas the term $R_{2,\varepsilon}$ can have arbitrary sign and scales with the diffusion coefficients $\lambda_{ij}$.

Hence, if the first term dominates, we can assure positivity of the integral. This follows with (4.10), (4.23) and choosing $C_{2}$ sufficiently large. Then Gronwall’s inequality yields a constant $K>0$ such that

[TABLE]

5 Conclusions

We have presented how to derive the system (2.22) in such a way that it automatically satisfies an entropy inequality and hence fulfills the second law of thermodynamics. There exist smooth solutions globally in time to this system if the smooth initial data are close enough to an equilibrium. In an asymptotic time regime we show the convergence to a parabolic limit system generalizing results on the single-component case.

Acknowledgments

The authors kindly acknowledge the financial support of this work by the Deutsche Forschungsgemeinschaft (DFG) in the frame of the International Research Training Group "Droplet Interaction Technologies" (DROPIT).

Appendix A Differential Operators and Matrix Algebra

We collect some definitions from vector analysis and matrix algebra which are frequently used in Sections 2–4.

For some vector $\mathbf{u}(\mathbf{x})=(u_{1}(\mathbf{x}),\dots,u_{n}(\mathbf{x}))^{\top}\in{\mathbb{R}}^{n}$ the (generalized) gradient a is defined as

[TABLE]

and for $\mathbf{v}(\mathbf{x})=(\mathbf{v}_{1}^{\top}(\mathbf{x}),\dots,\mathbf{v}_{n}^{\top}(\mathbf{x}))^{\top}\in{\mathbb{R}}^{nd}$ the (generalized) divergence is given by

[TABLE]

By $\otimes$ we denote the Kronecker product of two matrices, i.e., with $A\in{\mathbb{R}}^{m\times n},B\in{\mathbb{R}}^{p\times q}$

[TABLE]

We conclude with the following rules for the generalized gradient defined in (A.1).

Lemma 7**.**

For smooth functions $\alpha\colon{\mathbb{R}}^{d}\to{\mathbb{R}},\mathbf{a}\colon{\mathbb{R}}^{d}\to{\mathbb{R}}^{n},$ $\mathbf{b}\colon{\mathbb{R}}^{d}\to{\mathbb{R}}^{n},$ and $\mathbf{c}\colon{\mathbb{R}}^{n}\to{\mathbb{R}}^{n}$ it holds

[TABLE]

Proof A.1**.**

We have

[TABLE]

and

[TABLE]

Appendix B Global Classical Well-posedness of IVPs for Hyperbolic Balance Laws

Let the state space $G\subset{\mathbb{R}}^{n(d+1)}$ be open and denote by $U:{\mathbb{R}}^{d}\times[0,T)\to G$ the function that satisfies for some $T>0$ the IVP for the nonlinear system of balance laws given by

[TABLE]

Here $S:G\to{\mathbb{R}}^{n(d+1)}$ and $F_{\alpha}:G\to{\mathbb{R}}^{n(d+1)},\ \alpha=1,\dots,d$ are smooth functions with

[TABLE]

From now on we assume that $U$ can be split according to $U=(\mathbf{r}^{\top},\mathbf{m}^{\top})^{\top},$ with $\mathbf{r}\in{\mathbb{R}}^{n},\mathbf{m}\in{\mathbb{R}}^{nd}$. The system (B.1) can then be written as

[TABLE]

Theorem 8** ([25]).**

Let $s\geq s_{0}+1=\lfloor d/2\rfloor+2$ be an integer and $\hat{U}\in G$ be a constant equilibrium state such that the following conditions hold:

The Jacobian $\operatorname{D}_{\mathbf{m}}\mathbf{s}(\hat{U})\in{\mathbb{R}}^{nd\times nd}$ is invertible. 2. 2.

There exists a strictly convex smooth entropy function $\eta:\mathcal{G}\to{\mathbb{R}}$, defined in a convex, compact neighborhood $\mathcal{G}\subset G$ of $\hat{U}$, such that $\operatorname{D}^{2}\eta(U)\operatorname{D}F_{\alpha}(U)$ is symmetric for all $U\in\mathcal{G}$ and all $\alpha=1,\dots,d$. 3. 3.

There is a positive constant $c_{\mathcal{G}}$ such that for all $U\in\mathcal{G}$,

[TABLE]

where $|\cdot|$ denotes the Euclidean norm of a vector. 4. 4.

*The kernel $\ker(\operatorname{D}S(\hat{U}))$ of the Jacobian $\operatorname{D}S(\hat{U})$ contains no eigenvector of the matrix *

$\sum_{\alpha=1}^{d}\omega_{\alpha}\operatorname{D}F_{\alpha}(\hat{U})$, for any $\omega=(\omega_{1},\dots,\omega_{d})\in\mathbb{S}^{d-1}$.

Then there is a constant $c_{1}>0$ such that for $U_{0}\in H^{s}({\mathbb{R}}^{d})$ with

[TABLE]

the system of balance laws (B.2) with $U_{0}$ as its initial value has a unique global solution $U=U(\mathbf{x},t)\in C([0,T);H^{s}({\mathbb{R}}^{d}))$ satisfying

[TABLE]

for any $T>0$ and some $c_{2}>0$.