Rigorous Derivation of the Oberbeck–Boussinesq Approximation Revealing Unexpected Term

Abstract

We consider a general compressible viscous and heat conducting fluid confined between two parallel plates and heated from the bottom. The time evolution of the fluid is described by the Navier–Stokes–Fourier system considered in the regime of low Mach and Froude numbers suitably interrelated. Surprisingly and differently to the case of Neumann boundary conditions for the temperature, the asymptotic limit is identified as the Oberbeck–Boussinesq system supplemented with non-local boundary conditions for the temperature deviation.

1 Introduction

The Rayleigh–Bénard problem concerns the motion of a compressible fluid confined between two parallel plates heated from the bottom. For the sake of simplicity, we suppose the motion is space-periodic with respect to the horizontal variable. Accordingly, the underlying spatial domain may be identified with the flat torus in the horizontal plane:

$$\begin{aligned} \Omega = {\mathbb {T}}^{d - 1} \times (0,1),\ {\mathbb {T}}^{d-1} =\left( [0,1] \Big |_{\{ 0,1 \} } \right) ^{d - 1},\ d=2,3. \end{aligned}$$

The state of the fluid at a given time t and a spatial position $x \in \Omega$ is characterized by the mass density $\varrho = \varrho (t,x)$, the absolute temperature $\vartheta = \vartheta (t,x)$, and the velocity $\textbf{u}= \textbf{u}(t,x)$. We consider the motion in the asymptotically incompressible and stratified regime. Accordingly, the time evolution of the flow is governed by the scaled Navier–Stokes–Fourier (NSF) system:

The viscous stress is given by Newton’s rheological law

$$\begin{aligned} {\mathbb {S}}(\vartheta , \nabla _{x}\textbf{u}) = \mu (\vartheta ) \left( \nabla _{x}\textbf{u}+ \nabla _{x}^t \textbf{u}- \frac{2}{d} \textrm{div}_{x}\textbf{u}{\mathbb {I}} \right) + \eta (\vartheta ) \textrm{div}_{x}\textbf{u}{\mathbb {I}}, \end{aligned}$$

(1.6)

while the internal energy flux satisfies Fourier’s law

$$\begin{aligned} \textbf{q}(\vartheta , \nabla _{x}\vartheta ) = - \kappa (\vartheta ) \nabla _{x}\vartheta . \end{aligned}$$

(1.7)

Here, the more conventional internal energy (heat) equation

$$\begin{aligned} \partial _t (\varrho e(\varrho , \vartheta )) + \textrm{div}_{x}(\varrho e(\varrho , \vartheta ) \textbf{u}) + \textrm{div}_{x}\textbf{q} (\vartheta , \nabla _{x}\vartheta ) =\varepsilon ^2 {\mathbb {S}} (\vartheta , \nabla _{x}\textbf{u}): \nabla _{x}\textbf{u}- p (\varrho , \vartheta ) \textrm{div}_{x}\textbf{u}\end{aligned}$$

(1.8)

is replaced by the entropy balance equation (1.3). We point out that the equations (1.3), (1.8) are equivalent in the framework of classical solutions as long as the thermodynamic functions $p(\varrho , \vartheta )$, $e(\varrho , \vartheta )$, $s(\varrho , \vartheta )$ are interrelated through Gibbs’ equation

$$\begin{aligned} \vartheta D s = D e + p D \left( \frac{1}{\varrho } \right) . \end{aligned}$$

(1.9)

The entropy balance equation, however, is more convenient for the weak formulation introduced in Sect. 2. In particular, the weak solutions based on the entropy formulation enjoy several useful properties, notably the so-called weak strong uniqueness principle.

Besides Gibbs’ equation, we impose the hypothesis of thermodynamic stability written in the form

$$\begin{aligned} \frac{\partial p(\varrho , \vartheta ) }{\partial \varrho }> 0,\ \frac{\partial e(\varrho , \vartheta ) }{\partial \vartheta }> 0 \ \text{ for } \text{ all }\ \varrho , \vartheta > 0. \end{aligned}$$

(1.10)

The scaling in (1.1)–(1.3) represented by a small parameter $\varepsilon > 0$ corresponds to the Mach number $\textrm{Ma} = \varepsilon$ (the ratio between a “characteristic” fluid’s velocity and the local speed of sound) and the Froude number $\textrm{Fr} = \varepsilon ^{\frac{1}{2}}$ (the ratio between inertial and gravitational forces of the fluid), see e.g. the survey paper by Klein et al. [17] and our recent work with different scaling [4]. Our goal is to identify the asymptotic limit of solutions for $\varepsilon \rightarrow 0$.

1.1 Asymptotic limit

The zero-th order terms in the asymptotic limit are determined by the stationary problem

$$\begin{aligned} \nabla _{x}p(\varrho , \vartheta ) = \varepsilon \varrho \nabla _{x}G. \end{aligned}$$

(1.11)

Accordingly, a solution of the Navier–Stokes–Fourier system can be written in the form

$$\begin{aligned} \varrho _{\varepsilon }= \overline{\varrho } + \varepsilon \mathfrak {R}_\varepsilon ,\ \vartheta _{\varepsilon }= \overline{\vartheta } + \varepsilon \mathfrak {T}_\varepsilon , \end{aligned}$$

where $\overline{\varrho }$, $\overline{\vartheta }$ are positive constants and $\mathfrak {R}_\varepsilon \rightarrow \mathfrak {R}$, $\mathfrak {T}_\varepsilon \rightarrow \mathfrak {T}$ satisfy in the asymptotic limit $\varepsilon \rightarrow 0$ the so called Boussinesq relation

$$\begin{aligned} \frac{\partial p(\overline{\varrho }, \overline{\vartheta } ) }{\partial \varrho } \nabla _{x}\mathfrak {R} +\frac{\partial p(\overline{\varrho }, \overline{\vartheta } ) }{\partial \vartheta } \nabla _{x}\mathfrak {T} = \overline{\varrho } \nabla _{x}G. \end{aligned}$$

(1.12)

Without loss of generality, we suppose

meaning, in particular, the total mass of the fluid is constant independent of $\varepsilon$.

Anticipating convergence of the temperature deviations

$$\begin{aligned} \mathfrak {T}_\varepsilon = \frac{\vartheta _{\varepsilon }- \overline{\vartheta }}{\varepsilon } \rightarrow \mathfrak {T} \end{aligned}$$

we impose the boundary condition

$$\begin{aligned} \vartheta _{\varepsilon }|_{\partial \Omega } = \vartheta _B = \overline{\vartheta } + \varepsilon \Theta _B, \ \Theta _B = \Theta _B (x), \end{aligned}$$

(1.13)

where the perturbation $\Theta _B$ is not necessarily positive.

1.2 Limit system

Formally, the limit problem is expected to be the Oberbeck–Boussinesq (OB) system:

$$\begin{aligned} \textrm{div}_{x}\textbf{U}&= 0, \nonumber \\ \overline{\varrho } \Big ( \partial _t \textbf{U}+ \textbf{U}\cdot \nabla _{x}\textbf{U}\Big ) + \nabla _{x}\Pi&= \textrm{div}_{x}{\mathbb {S}}(\overline{\vartheta }, \nabla _{x}\textbf{U}) +\mathfrak {R} \nabla _{x}G, \nonumber \\ \overline{\varrho } c_p(\overline{\varrho }, \overline{\vartheta } ) \left( \partial _t \mathfrak {T} + \textbf{U}\cdot \nabla _{x}\mathfrak {T} \right) -\overline{\varrho } \ \overline{\vartheta } \alpha (\overline{\varrho }, \overline{\vartheta } ) \textbf{U}\cdot \nabla _{x}G&= \kappa (\overline{\vartheta }) \Delta _x\mathfrak {T}, \end{aligned}$$

(1.14)

where

$$\begin{aligned} \alpha (\overline{\varrho }, \overline{\vartheta } )&\equiv \frac{1}{\overline{\varrho }} \frac{\partial p(\overline{\varrho }, \overline{\vartheta } ) }{\partial \vartheta } \left( \frac{\partial p(\overline{\varrho }, \overline{\vartheta } ) }{\partial \varrho } \right) ^{-1},\nonumber \\ c_p (\overline{\varrho }, \overline{\vartheta } )&\equiv \frac{\partial e(\overline{\varrho }, \overline{\vartheta } ) }{\partial \vartheta } + \overline{\varrho }^{-1} \overline{\vartheta } \alpha (\overline{\varrho }, \overline{\vartheta } ) \frac{\partial p(\overline{\varrho }, \overline{\vartheta } ) }{\partial \vartheta } \end{aligned}$$

(1.15)

stand for the coefficient of thermal expansion, and the specific heat at constant pressure, respectively. The quantities $\mathfrak {R}$, $\mathfrak {T}$ satisfy the Boussinesq relation (1.12).

We refer to the survey by Zeytounian [21] and the list of references therein for a formal derivation of the Oberbeck–Boussinesq system. A rigorous proof of convergence of a family of (weak) solutions $(\varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon )_{\varepsilon > 0}$ to the Navier–Stokes–Fourier system, meaning

$$\begin{aligned} \textbf{u}_\varepsilon \rightarrow \textbf{U}, \ \frac{\vartheta _{\varepsilon }- \overline{\vartheta }}{\varepsilon } \rightarrow \mathfrak {T} \ \text{ in } \text{ some } \text{ sense, } \end{aligned}$$

was obtained in [11] (see also [12, Chapter 5]) for the conservative boundary conditions

$$\begin{aligned} \textbf{u}_\varepsilon \cdot \textbf{n}|_{\partial \Omega }&= 0,\ ({\mathbb {S}} (\vartheta _{\varepsilon }, \nabla _{x}\textbf{u}_\varepsilon ) \cdot \textbf{n}) \times \textbf{n}|_{\partial \Omega } = 0, \nonumber \\ \nabla _{x}\vartheta _{\varepsilon }\cdot \textbf{n}|_{\partial \Omega }&= 0. \end{aligned}$$

(1.16)

Note that the conservative boundary conditions (1.16) imply that the average

is a constant of motion. In particular, the Boussinesq relation (1.12) can be written as

$$\begin{aligned} \frac{\partial p(\overline{\varrho }, \overline{\vartheta } ) }{\partial \varrho } \mathfrak {R} +\frac{\partial p(\overline{\varrho }, \overline{\vartheta } ) }{\partial \vartheta } \mathfrak {T} = \overline{\varrho } G \end{aligned}$$

(1.17)

recalling the convention $\int _{\Omega } {G } \ \,\textrm{d} {x} = 0$. The form (1.17) of the Boussinesq relation is rather misleadingly interpreted in the sense that the Oberbeck–Boussinesq approximation can be derived from the Navier–Stokes–Fourier system on condition that the density is supposed to be an affine function of the temperature, cf. the discussion in Maruyama [18].

The present paper can be seen as the first attempt to address the same singular limit problem under the inhomogeneous Dirichlet boundary conditions imposed on the temperature, physically relevant for the well-known Rayleigh–Bénard convection problem. The result is stated in the framework of the general theory of weak solutions to the primitive Navier–Stokes–Fourier system with inhomogeneous boundary conditions developed recently in [8] and the monograph [13].

Rather surprisingly, we show that the limit system is in fact different if the Dirichlet boundary conditions (1.4), (1.5) are imposed. In accordance with (1.4), (1.13), we consider the scaled Navier–Stokes–Fourier system with the boundary conditions

$$\begin{aligned} \textbf{u}_\varepsilon |_{\partial \Omega }= 0,\ \vartheta _{\varepsilon }|_{\partial \Omega } = \vartheta _{B}= \overline{\vartheta } + \varepsilon \Theta _B. \end{aligned}$$

(1.18)

Imposing (1.18) in place of (1.16) drastically changes the behaviour of the fluid as the resulting system is no longer energetically closed and the dynamics is driven by the temperature gradient. Indeed the total mass conservation discussed by several authors (Barletta et al. [2, 3], Maruyama [18]) must be enforced through the condition

$$\begin{aligned} \int _{\Omega } { \mathfrak {R}_\varepsilon } \ \,\textrm{d} {x} = 0 \ \Rightarrow \ \int _{\Omega } { \mathfrak {R} } \ \,\textrm{d} {x} = 0 \end{aligned}$$

in (1.12). Accordingly, the conventional Boussinesq relation (1.17) must be replaced by

$$\begin{aligned} \frac{\partial p(\overline{\varrho }, \overline{\vartheta } ) }{\partial \varrho } \mathfrak {R} +\frac{\partial p(\overline{\varrho }, \overline{\vartheta } ) }{\partial \vartheta } \mathfrak {T} = \overline{\varrho } G + \chi (t), \end{aligned}$$

(1.19)

where, in accordance with our convention $\int _{\Omega } {G} \ \,\textrm{d} {x} = 0$,

(1.20)

Relation (1.20) indicates that the limit problem could involve “non-local” terms. Indeed we show that a family of (weak) solutions to the NSF system $(\varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon )_{\varepsilon > 0}$ admits a limit

$$\begin{aligned} \textbf{u}_\varepsilon \rightarrow \textbf{U}, \ \frac{\vartheta _{\varepsilon }- \overline{\vartheta }}{\varepsilon } \rightarrow \mathfrak {T},\ \frac{\varrho _{\varepsilon }- \overline{\varrho }}{\varepsilon } \rightarrow \mathfrak {R} \end{aligned}$$

solving a modified Oberbeck–Bousinesq system

(1.21)

together with the Boussinesq relation

$$\begin{aligned} \frac{\partial p(\overline{\varrho }, \overline{\vartheta } ) }{\partial \varrho } \nabla _{x}\mathfrak {R} +\frac{\partial p(\overline{\varrho }, \overline{\vartheta } ) }{\partial \vartheta } \nabla _{x}\mathfrak {T} = \overline{\varrho } \nabla _{x}G,\ \int _{\Omega } { \mathfrak {R} } \ \,\textrm{d} {x} = 0, \end{aligned}$$

(1.22)

and the boundary conditions

$$\begin{aligned} \textbf{U}|_{\partial \Omega } = 0,\ \mathfrak {T}|_{\partial \Omega } = \Theta _B. \end{aligned}$$

(1.23)

To simplify even more, using relation (1.19) to replace $\mathfrak {R}$ with $\mathfrak {T}$ in the momentum equation while observing that two gradient terms hide into $\Pi$, the system (1.21)–(1.23) turns into

(1.24)

Finally, performing a simple change of variables

we may rewrite system (1.21)–(1.23) in the form of conventional Oberbeck–Boussinesq system

Here and hereafter, we always suppose that the coefficient of thermal expansion $\alpha (\overline{\varrho }, \overline{\vartheta })$ at the reference state $(\overline{\varrho }, \overline{\vartheta })$ is strictly positive. In particular, it follows directly from (1.15) and (1.10) that $\lambda \in (0,1)$. Parabolic equations with the non-local boundary terms similar to (1.27) have been recognized by Day [10] in the context of some models in thermoelasticity, and subsequently studied by a number of authors, see Chen and Liu [9], Friedman [14], Gladkov and Nikitin [16], Pao [20], among others.

We point out that our results justify the Oberbeck–Boussinesq approximation as a rigorous description of the limit behaviour of solutions to the complete Navier–Stokes–Fourier system only in the regime where the temperature deviation from the constant state $\overline{\vartheta }$ is small of order $\varepsilon$. Indeed the relevance of the Oberbeck–Boussinesq approximation for the Rayleigh–Bénard convection with large deviation of the boundary temperature has been questioned by several authors, see Bormann [7], Fröhlich, Laure, and Peyret [15], Nadolin [19], among others.

1.3 The strategy of the convergence proof

Our approach is based on the concept of weak solutions to the NSF system with energetically open boundary conditions developed recently in [8, 13]. In particular, the relative energy inequality based on the ballistic energy balance is used to measure the distance between the solutions of the primitive and target systems. In contrast with [11, 12, Chapter 5], strong convergence is obtained with certain explicit estimates of the error depending on how close are the initial data to their limit values.

The paper is organized as follows. We start with the concept of weak solutions for the NSF system with Dirichlet boundary conditions introduced in [8]. In particular, we recall the ballistic energy and the associated relative energy inequality in Sect. 2. Then, in Sect. 3, we record the available results on solvability of the OB system in the framework of strong solutions. The main results on convergence to the target OB system are stated in Sect. 4. The rest of the paper is then devoted to the proof of convergence. In Sect. 5, we derive the basic energy estimates that control the amplitude of the fluid velocity as well as the distance of the density and the temperature profiles from their limit values independent of the scaling parameter $\varepsilon$. The proof of convergence to the OB system is completed in Sect. 6 by means of the relative energy inequality.

2 Weak Solutions to the Primitive System

Following [8, 13], we introduce the concept of weak solutions to the NSF system.

Definition 2.1

(Weak solution to the NSF system). We say that a trio $(\varrho , \vartheta , \textbf{u})$ is a weak solution of the NSF system (1.1)–(1.7), with the initial data

$$\begin{aligned} \varrho (0, \cdot ) = \varrho _0,\ (\varrho \textbf{u}) (0, \cdot ) = \varrho _0 \textbf{u}_0,\ (\varrho s)(0, \cdot ) = \varrho _0 s(\varrho _0, \vartheta _0), \end{aligned}$$

if the following holds:

• The solution belongs to the regularity class:

  $$\begin{aligned} \varrho&\in L^\infty (0,T; L^\gamma (\Omega )) \ \text{ for } \text{ some }\ \gamma> 1,\ \varrho \ge 0 \ \text{ a.a. } \text{ in }\ (0,T) \times \Omega , \nonumber \\ \textbf{u}&\in L^2(0,T; W^{1,2}_0 (\Omega ; {\mathbb {R}}^d)), \nonumber \\ \vartheta ^{\beta /2} ,\ \log (\vartheta )&\in L^2(0,T; W^{1,2}(\Omega )) \ \text{ for } \text{ some }\ \beta \ge 2,\ \vartheta > 0 \ \text{ a.a. } \text{ in }\ (0,T) \times \Omega , \nonumber \\ (\vartheta - \vartheta _{B})&\in L^2(0,T; W^{1,2}_0 (\Omega )). \end{aligned}$$

  (2.1)

• The equation of continuity (1.1) is satisfied in the sense of distributions, more specifically,

  $$\begin{aligned} \int _0^T \int _{\Omega } { \Big [ \varrho \partial _t \varphi + \varrho \textbf{u}\cdot \nabla _{x}\varphi \Big ] } \ \,\textrm{d} {x} \,\textrm{d} t&= - \int _{\Omega } { \varrho (0) \varphi (0, \cdot ) } \ \,\textrm{d} {x} \end{aligned}$$

  (2.2)

  for any $\varphi \in C^1_c([0,T) \times \overline{\Omega } )$.

• The momentum equation (1.2) is satisfied in the sense of distributions,

  $$\begin{aligned}&\int _0^T \int _{\Omega } { \left[ \varrho \textbf{u}\cdot \partial _t \varvec{\varphi }+ \varrho \textbf{u}\otimes \textbf{u}: \nabla _{x}\varvec{\varphi }+ \frac{1}{\varepsilon ^2} p(\varrho , \vartheta ) \textrm{div}_{x}\varvec{\varphi }\right] } \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\quad = \int _0^T \int _{\Omega } { \left[ {\mathbb {S}} (\vartheta , \nabla _{x}\textbf{u}) : \nabla _{x}\varvec{\varphi }- \frac{1}{\varepsilon } \varrho \nabla _{x}G \cdot \varvec{\varphi }\right] } \ \,\textrm{d} {x} \,\textrm{d} t- \int _{\Omega } { \varrho _0 \textbf{u}_0 \cdot \varvec{\varphi }(0, \cdot ) } \ \,\textrm{d} {x} \end{aligned}$$

  (2.3)

  for any $\varvec{\varphi }\in C^1_c([0, T) \times \Omega ; {\mathbb {R}}^d)$.

• The entropy balance (1.3) is replaced by the inequality

  $$\begin{aligned}&- \int _0^T \int _{\Omega } { \left[ \varrho s(\varrho , \vartheta ) \partial _t \varphi + \varrho s (\varrho ,\vartheta ) \textbf{u}\cdot \nabla _{x}\varphi + \frac{\textbf{q} (\vartheta , \nabla _{x}\vartheta )}{\vartheta } \cdot \nabla _{x}\varphi \right] } \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\quad \ge \int _0^T \int _{\Omega } { \frac{\varphi }{\vartheta } \left[ \varepsilon ^2 {\mathbb {S}}(\vartheta , \nabla _{x}\textbf{u}) : \nabla _{x}\textbf{u}- \frac{\textbf{q} (\vartheta , \nabla _{x}\vartheta ) \cdot \nabla _{x}\vartheta }{\vartheta } \right] } \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\qquad + \int _{\Omega } { \varrho _0 s(\varrho _0, \vartheta _0) \varphi (0, \cdot ) } \ \,\textrm{d} {x} \end{aligned}$$

  (2.4)

  for any $\varphi \in C^1_c([0, T) \times \Omega )$, $\varphi \ge 0$.

• The ballistic energy balance

  $$\begin{aligned}&- \int _0^T \partial _t \psi \int _{\Omega } { \left[ \varepsilon ^2 \frac{1}{2} \varrho |\textbf{u}|^2 + \varrho e(\varrho , \vartheta ) - \tilde{\vartheta }\varrho s(\varrho , \vartheta ) \right] } \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\qquad + \int _0^T \psi \int _{\Omega } { \frac{\tilde{\vartheta }}{\vartheta } \left[ \varepsilon ^2 {\mathbb {S}}(\vartheta , \nabla _{x}\textbf{u}): \nabla _{x}\textbf{u}-\frac{\textbf{q}(\vartheta , \nabla _{x}\vartheta ) \cdot \nabla _{x}\vartheta }{\vartheta } \right] } \ \,\textrm{d} {x}\,\textrm{d} t\nonumber \\&\quad \le \int _0^T \psi \int _{\Omega } { \left[ \varepsilon \varrho \textbf{u}\cdot \nabla _{x}G - \varrho s(\varrho , \vartheta ) \partial _t \tilde{\vartheta }- \varrho s(\varrho , \vartheta ) \textbf{u}\cdot \nabla _{x}\tilde{\vartheta }- \frac{\textbf{q}(\vartheta , \nabla _{x}\vartheta )}{\vartheta } \cdot \nabla _{x}\tilde{\vartheta }\right] } \ \,\textrm{d} {x}\nonumber \\&\qquad + \psi (0) \int _{\Omega } { \left[ \frac{1}{2} \varepsilon ^2 \varrho _0 |\textbf{u}_0|^2 + \varrho _0 e(\varrho _0, \vartheta _0) - \tilde{\vartheta }(0, \cdot ) \varrho _0 s(\varrho _0, \vartheta _0) \right] } \ \,\textrm{d} {x} \end{aligned}$$

  (2.5)

  holds for any $\psi \in C^1_c[0; T)$, $\psi \ge 0$, and any $\tilde{\vartheta }\in C^1([0; T) \times \overline{\Omega })$ satisfying

  $$\begin{aligned} \tilde{\vartheta }> 0,\ \tilde{\vartheta }|_{\partial \Omega } = \vartheta _{B}. \end{aligned}$$

Although quite general, the weak solutions in the sense of Definition 2.1 comply with the weak–strong uniqueness principle, meaning they coincide with the strong solution as long as the latter exists, see [13, Chapter 4].

2.1 Relative energy inequality

Following [8, 13], we introduce the scaled relative energy

$$\begin{aligned} \begin{aligned}&E_\varepsilon \left( \varrho , \vartheta , {\textbf {u}}\Big | \tilde{\varrho }, \tilde{\vartheta }, {\tilde{{\textbf {u}}}}\right) = \frac{1}{2}\varrho |{\textbf {u}}- {\tilde{{\textbf {u}}}}|^2 + \\&\qquad \frac{1}{\varepsilon ^2} \left[ \varrho e - \tilde{\vartheta }\Big (\varrho s - \tilde{\varrho }s(\tilde{\varrho }, \tilde{\vartheta }) \Big ) -\Big ( e(\tilde{\varrho }, \tilde{\vartheta }) - \tilde{\vartheta }s(\tilde{\varrho }, \tilde{\vartheta }) +\frac{p(\tilde{\varrho }, \tilde{\vartheta })}{\tilde{\varrho }} \Big ) (\varrho - \tilde{\varrho }) - \tilde{\varrho }e (\tilde{\varrho }, \tilde{\vartheta }) \right] . \end{aligned} \end{aligned}$$

The hypothesis of thermodynamic stability (1.10) can be equivalently rephrased as (strict) convexity of the total energy expressed with respect to the conservative entropy variables

$$\begin{aligned} E_\varepsilon \Big ( \varrho , S = \varrho s(\varrho , \vartheta ), \textbf{m}= \varrho \textbf{u}\Big ) \equiv \frac{1}{2} \frac{|\textbf{m}|^2}{\varrho } + \frac{1}{\varepsilon ^2} \varrho e(\varrho , S), \end{aligned}$$

whereas the relative energy can be written as Bregmann distance

$$\begin{aligned} E_\varepsilon \left( \varrho , S, \textbf{m}\Big | \tilde{\varrho }, \widetilde{S}, \widetilde{\textbf{m}}\right)&= E_\varepsilon (\varrho , S, \textbf{m}) - \left\langle \partial _{\varrho , S, \textbf{m}} E_\varepsilon (\tilde{\varrho }, \widetilde{S}, \widetilde{\textbf{m}}) ; (\varrho - \tilde{\varrho }, S - \widetilde{S}, \textbf{m}- \widetilde{\textbf{m}}) \right\rangle \\&\quad - E_\varepsilon (\tilde{\varrho }, \widetilde{S}, \widetilde{\textbf{m}}) \end{aligned}$$

associated to the convex functional $(\varrho , S, \textbf{m}) \mapsto E_\varepsilon (\varrho , S, \textbf{m})$, see [13, Chapter 3, Section 3.1]. Finally, as stated in [8], any weak solution in the sense of Definition 2.1 satisfies the relative energy inequality

$$\begin{aligned}&\left[ \int _{\Omega } { E_\varepsilon \left( \varrho , \vartheta , {\textbf {u}}\Big | \tilde{\varrho }, \tilde{\vartheta }, {\tilde{{\textbf {u}}}}\right) } \ \,\text {d} {x} \right] _{t = 0}^{t = \tau } \nonumber \\&\qquad + \int _0^\tau \int _{\Omega } { \frac{\tilde{\vartheta }}{\vartheta } \left( {\mathbb {S}} (\vartheta , {\mathbb {D}}_x{\textbf {u}}) : {\mathbb {D}}_x{\textbf {u}}+ \frac{1}{\varepsilon ^2} \frac{\kappa (\vartheta ) |\nabla _{x}\vartheta |^2 }{\vartheta } \right) } \ \,\text {d} {x} \,\text {d} t\nonumber \\ {}&\quad \le - \frac{1}{\varepsilon ^2} \int _0^\tau \int _{\Omega } { \left( \varrho (s - s(\tilde{\varrho }, \tilde{\vartheta })) \partial _t \tilde{\vartheta }+ \varrho (s - s(\tilde{\varrho }, \tilde{\vartheta })) {\textbf {u}}\cdot \nabla _{x}\tilde{\vartheta }- \left( \frac{\kappa (\vartheta ) \nabla _{x}\vartheta }{\vartheta } \right) \cdot \nabla _{x}\tilde{\vartheta }\right) } \ \,\text {d} {x} \,\text {d} t\nonumber \\ {}&\qquad - \int _0^\tau \int _{\Omega } { \Big [ \varrho ({\textbf {u}}- {\tilde{{\textbf {u}}}}) \otimes ({\textbf {u}}- {\tilde{{\textbf {u}}}}) + \frac{1}{\varepsilon ^2} p(\varrho , \vartheta ) {\mathbb {I}} - {\mathbb {S}}(\vartheta , {\mathbb {D}}_x{\textbf {u}}) \Big ] : {\mathbb {D}}_x{\tilde{{\textbf {u}}}}} \ \,\text {d} {x} \,\text {d} t\nonumber \\ {}&\qquad + \int _0^\tau \int _{\Omega } { \varrho \left[ \frac{1}{\varepsilon } \nabla _{x}G - \partial _t {\tilde{{\textbf {u}}}}- ({\tilde{{\textbf {u}}}}\cdot \nabla _{x}) {\tilde{{\textbf {u}}}}\right] \cdot ({\textbf {u}}- {\tilde{{\textbf {u}}}}) } \ \,\text {d} {x} \,\text {d} t\nonumber \\ {}&\qquad + \frac{1}{\varepsilon ^2} \int _0^\tau \int _{\Omega } { \left[ \left( 1 - \frac{\varrho }{\tilde{\varrho }} \right) \partial _t p(\tilde{\varrho }, \tilde{\vartheta }) - \frac{\varrho }{\tilde{\varrho }} {\textbf {u}}\cdot \nabla _{x}p(\tilde{\varrho }, \tilde{\vartheta }) \right] } \ \,\text {d} {x} \,\text {d} t \end{aligned}$$

(2.6)

for a.a. $\tau > 0$ and any trio of continuously differentiable functions $(\tilde{\varrho }, \tilde{\vartheta }, {\tilde{\textbf{u}}})$ satisfying

$$\begin{aligned} \tilde{\varrho }> 0,\ \tilde{\vartheta }> 0,\ \tilde{\vartheta }|_{\partial \Omega } = \vartheta _{B}, \ {\tilde{\textbf{u}}}|_{\partial \Omega } = 0. \end{aligned}$$

(2.7)

Here, we have denoted ${\mathbb {D}}_x\textbf{u}=\frac{1}{2}(\nabla _{x}\textbf{u}+\nabla _{x}^t \textbf{u})$.

2.2 Constitutive relations

The existence theory developed in [8, 13] is conditioned by certain restrictions imposed on the constitutive relations (state equations) similar to those introduced in the monograph [12, Chapters 1,2]. Specifically, the equation of state takes the form

$$\begin{aligned} p(\varrho , \vartheta ) = p_{\textrm{m}} (\varrho , \vartheta ) + p_{\textrm{rad}}(\vartheta ), \end{aligned}$$

(2.8)

where $p_{\textrm{m}}$ is the pressure of a general monoatomic gas,

$$\begin{aligned} p_{\textrm{m}} (\varrho , \vartheta ) = \frac{2}{3} \varrho e_{\textrm{m}}(\varrho , \vartheta ), \end{aligned}$$

(2.9)

enhanced by the radiation pressure

$$\begin{aligned} p_{\textrm{rad}}(\vartheta ) = \frac{a}{3} \vartheta ^4,\ a > 0. \end{aligned}$$

Accordingly, the internal energy reads

$$\begin{aligned} e(\varrho , \vartheta ) = e_{\textrm{m}}(\varrho , \vartheta ) + e_{\textrm{rad}}(\varrho , \vartheta ), \ e_{\textrm{rad}}(\varrho , \vartheta ) = \frac{a}{\varrho } \vartheta ^4. \end{aligned}$$

Moreover, using several physical principles it was shown in [12, Chapter 1]:

• Gibbs’ relation together with (2.9) yield

  $$\begin{aligned} p_{\textrm{m}} (\varrho , \vartheta ) = \vartheta ^{\frac{5}{2}} P \left( \frac{\varrho }{\vartheta ^{\frac{3}{2}} } \right) \end{aligned}$$

  for a certain $P \in C^1[0,\infty )$. Consequently,

  $$\begin{aligned} p(\varrho , \vartheta ) = \vartheta ^{\frac{5}{2}} P \left( \frac{\varrho }{\vartheta ^{\frac{3}{2}}} \right) + \frac{a}{3} \vartheta ^4,\ e(\varrho , \vartheta ) = \frac{3}{2} \frac{\vartheta ^{\frac{5}{2}} }{\varrho } P \left( \frac{\varrho }{\vartheta ^{\frac{3}{2}}} \right) + \frac{a}{\varrho } \vartheta ^4, \ a > 0. \end{aligned}$$

  (2.10)

• Hypothesis of thermodynamic stability (1.10) expressed in terms of P gives rise to

  $$\begin{aligned} P(0) = 0,\ P'(Z)> 0 \ \text{ for }\ Z \ge 0, \ 0 < \frac{ \frac{5}{3} P(Z) - P'(Z) Z }{Z} \le C \ \text{ for }\ Z > 0. \end{aligned}$$

  (2.11)

  In particular, the function $Z \mapsto P(Z)/ Z^{\frac{5}{3}}$ is decreasing, and we suppose

  $$\begin{aligned} \lim _{Z \rightarrow \infty } \frac{ P(Z) }{Z^{\frac{5}{3}}} = p_\infty > 0. \end{aligned}$$

  (2.12)

• Accordingly, the associated entropy takes the form

  $$\begin{aligned} s(\varrho , \vartheta ) = s_{\textrm{m}}(\varrho , \vartheta ) + s_{\textrm{rad}}(\varrho , \vartheta ), \ s_{\textrm{m}} (\varrho , \vartheta ) = {\mathcal {S}} \left( \frac{\varrho }{\vartheta ^{\frac{3}{2}} } \right) , \ s_{\textrm{rad}}(\varrho , \vartheta ) = \frac{4a}{3} \frac{\vartheta ^3}{\varrho }, \end{aligned}$$

  (2.13)

  where

  $$\begin{aligned} {\mathcal {S}}'(Z) = -\frac{3}{2} \frac{ \frac{5}{3} P(Z) - P'(Z) Z }{Z^2} < 0. \end{aligned}$$

  (2.14)

• In addition, we impose the Third law of thermodynamics, cf. Belgiorno [5, 6], requiring the entropy to vanish when the absolute temperature approaches zero,

  $$\begin{aligned} \lim _{Z \rightarrow \infty } {\mathcal {S}}(Z) = 0. \end{aligned}$$

  (2.15)

Finally, we suppose the transport coefficients are continuously differentiable functions satisfying

$$\begin{aligned} 0< \underline{\mu }(1 + \vartheta )&\le \mu (\vartheta ),\ |\mu '(\vartheta )| \le \overline{\mu }, \nonumber \\ 0 \le \eta (\vartheta )&\le \overline{\eta }(1 + \vartheta ), \nonumber \\ 0 < \underline{\kappa } (1 + \vartheta ^\beta )&\le \kappa (\vartheta ) \le \overline{\kappa }(1 + \vartheta ^\beta ), \ \text{ where }\ \beta > 6. \end{aligned}$$

(2.16)

As a consequence of the above hypotheses, we get the following estimates:

$$\begin{aligned} \varrho ^{\frac{5}{3}} + \vartheta ^4 \lesssim \varrho e(\varrho , \vartheta )&\lesssim 1+ \varrho ^{\frac{5}{3}} + \vartheta ^4, \end{aligned}$$

(2.17)

$$\begin{aligned} s_{\textrm{m}}(\varrho , \vartheta )&\lesssim \left( 1 + |\log (\varrho )| + [\log (\vartheta )]^+ \right) , \end{aligned}$$

(2.18)

see [12, Chapter 3, Section 3.2]. Here and hereafter, the symbol $a \lesssim b$ means there is a positive constant $C> 0$ such that $a \le C b$.

We report the existence result proved in [8, Theorem 4.2].

Proposition 2.2

(Primitive system, global existence). Let the thermodynamic functions p, e, s and the transport coefficients $\mu$, $\eta$, $\kappa$ satisfy the hypotheses (2.8)–(2.16). Let $\vartheta _B \in C^2(\partial \Omega )$. Suppose the initial data $\varrho _0$, $\vartheta _0$, $\textbf{u}_0$ satisfy

$$\begin{aligned} \varrho _0> 0,\ \vartheta _0 > 0,\ \int _{\Omega } { E_\varepsilon (\varrho _0, \vartheta _0, \textbf{u}_0 ) } \ \,\textrm{d} {x} < \infty . \end{aligned}$$

Then for any $T > 0$, the Navier–Stokes–Fourier system (1.1)–(1.7) admits a weak solution $(\varrho , \vartheta , \textbf{u})$ in $(0,T) \times \Omega$ in the sense specified in Definition 2.1.

Remark 2.3

The constitutive restrictions imposed through (2.8)–(2.16) are necessary for the global in time existence result stated in Proposition 2.2. They are of technical nature, and, in our opinion, without significant impact of the asymptotic limit stated below. In particular, the relations (2.8), (2.9) can be replaced by more general ones specified in the monograph [12, Chapter 1, Sections 1.4.2, 1.4.3]. The radiation pressure component is essential only for the existence theory eliminating hypothetical temperature oscillations in the vacuum zones. The fact that its presence has no influence on the form of the asymptotic problem can be demonstrated by considering an extra scaling $a = a(\varepsilon ) \rightarrow 0$ as $\varepsilon \rightarrow 0$.

3 Strong Solutions to the Target System

Our analysis requires the existence of regular solutions to the Oberbeck–Boussinesq system (1.25)–(1.27). The relevant result was proved in [1, Theorem 2.3, Theorem 3.1].

Proposition 3.1

(Strong solutions to target system) Suppose that

$$\begin{aligned} G \in W^{1, \infty }({\Omega }),\ \Theta _B \in C^2(\overline{\Omega }), \end{aligned}$$

(3.1)

and

(3.2)

Then there exists $T_{\textrm{max}} > 0$, $T_{\textrm{max}} = \infty$ if $d = 2$, such that the OB system (1.25)–(1.27) with the initial data

$$\begin{aligned} \textbf{U}(0, \cdot ) = \textbf{U}_0,\ \Theta (0, \cdot ) = \Theta _0, \end{aligned}$$

admits a strong solution $\textbf{U}$, $\Theta$ in the regularity class

$$\begin{aligned}&\textbf{U}\in L^p(0,T; W^{2,p}(\Omega ; {\mathbb {R}}^d)), \ \partial _t \textbf{U}\in L^p(0,T; L^{p}(\Omega ; {\mathbb {R}}^d)), \ \Pi \in L^p(0,T; W^{1,p}(\Omega )), \nonumber \\&\Theta \in L^p(0,T; W^{2,p}(\Omega )), \ \partial _t \Theta \in L^p(0,T; L^{p}(\Omega ; {\mathbb {R}}^d)) \end{aligned}$$

(3.3)

for any $1 \le p < \infty$ and any $0< T < T_{\textrm{max}}$.

Remark 3.2

Strictly speaking, the existence result in [1] requires, in addition to (3.1),

$$\begin{aligned} \Delta _xG = 0. \end{aligned}$$

(3.4)

On the one hand, in applications, G represents the gravitational potential therefore (3.4) is automatically satisfied. On the other hand, the proof presented in [1] can be easily modified to accommodate the general case.

Given the parabolic character of the OB system, the strong solutions are in fact classical if higher regularity of the data is required, cf. [1, Theorem 4.1].

4 Asymptotic Limit, Main Result

We are ready to state our main result concerning the singular limit $\varepsilon \rightarrow 0$ in the primitive NSF system.

Note that hypotheses (4.3)–(4.5) correspond to well prepared data. The rest of the paper is devoted to the proof of Theorem 4.1.

5 Basic Energy Estimates

In order to perform the limit claimed in Theorem 4.1, we need bounds on the sequence $(\varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon )_{\varepsilon > 0}$ independent of $\varepsilon$. We start by introducing the notation borrowed from [12] distinguishing the “essential” and “residual” range of the thermostatic variables $(\varrho , \vartheta )$. Specifically, given a compact set

$$\begin{aligned} K \subset \left\{ (\varrho , \vartheta ) \in {\mathbb {R}}^2 \ \Big | \ \varrho> 0, \vartheta > 0 \right\} \end{aligned}$$

and $\varepsilon > 0$, we introduce

$$\begin{aligned} g_{\textrm{ess}} = g \mathbb {1}_{(\varrho _{\varepsilon }, \vartheta _{\varepsilon }) \in K}, \ g_{\textrm{res}} = g - g_{\textrm{ess}} = g \mathbb {1}_{(\varrho _{\varepsilon }, \vartheta _{\varepsilon }) \in {\mathbb {R}}^2 \setminus K} \end{aligned}$$

for any measurable $g = g(t,x)$. Note carefully that this decomposition depends on $\varepsilon$. As a matter of fact, the characteristic function $\mathbb {1}_{(\varrho _{\varepsilon }, \vartheta _{\varepsilon }) \in K}$ can be replaced by its smooth regularization by a suitable convolution kernel.

In the subsequent analysis, we consider

$$\begin{aligned} K = \overline{{\mathcal {U}}(\overline{\varrho }, \overline{\vartheta })} \subset (0, \infty )^2,\ {{\mathcal {U}} (\overline{\varrho }, \overline{\vartheta } )} \ \text{- } \text{ an } \text{ open } \text{ neighborhood } \text{ of }\ (\overline{\varrho }, \overline{\vartheta }). \end{aligned}$$

As shown in [12, Chapter 5, Lemma 5.1], the relative energy enjoys the following coercivity properties:

$$\begin{aligned} E_{\varepsilon } \left( \varrho , \vartheta , \textbf{u}\Big | \tilde{\varrho }, \tilde{\vartheta }, {\tilde{\textbf{u}}}\right) \ge C \left( \frac{ |\varrho - \tilde{\varrho }|^2 }{\varepsilon ^2} + \frac{ |\vartheta - \tilde{\vartheta }|^2 }{\varepsilon ^2} + |\textbf{u}- {\tilde{\textbf{u}}}|^2 \right) \end{aligned}$$

if $(\varrho , \vartheta ) \in K = \overline{{\mathcal {U}}(\overline{\varrho }, \overline{\vartheta })}$, $(\tilde{\varrho }, \tilde{\vartheta }) \in {\mathcal {U}}(\overline{\varrho }, \overline{\vartheta })$,

$$\begin{aligned} E_{\varepsilon } \left( \varrho , \vartheta , \textbf{u}\Big | \tilde{\varrho }, \tilde{\vartheta }, {\tilde{\textbf{u}}}\right) \ge C \left( \frac{1}{\varepsilon ^2} + \frac{1}{\varepsilon ^2} \varrho e(\varrho , \vartheta ) + \frac{1}{\varepsilon ^2} \varrho |s(\varrho , \vartheta )| + \varrho |\textbf{u}|^2 \right) \end{aligned}$$

whenever $(\varrho , \vartheta ) \in R^2 {\setminus } \overline{{\mathcal {U}}(\overline{\varrho }, \overline{\vartheta })}$, $(\tilde{\varrho }, \tilde{\vartheta }) \in {\mathcal {U}}(\overline{\varrho }, \overline{\vartheta })$. The constant C depends on K and the distance

$$\begin{aligned} \sup _{t,x} \textrm{dist} \left[ (\tilde{\varrho }(t,x), \tilde{\vartheta }(t,x) ); \partial K \right] . \end{aligned}$$

In other words,

$$\begin{aligned} E_{\varepsilon } \left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon \Big | \tilde{\varrho }, \tilde{\vartheta }, {\tilde{\textbf{u}}}\right) _{\textrm{ess}}&\ge C \left( \frac{ |\varrho _{\varepsilon }- \tilde{\varrho }|^2 }{\varepsilon ^2} +\frac{ |\vartheta _{\varepsilon }- \tilde{\vartheta }|^2 }{\varepsilon ^2} + |\textbf{u}_\varepsilon - {\tilde{\textbf{u}}}|^2 \right) _{\textrm{ess}}, \end{aligned}$$

(5.1)

$$\begin{aligned} E_{\varepsilon } \left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon \Big | \tilde{\varrho }, \tilde{\vartheta }, {\tilde{\textbf{u}}}\right) _{\textrm{res}}&\ge C \left( \frac{1}{\varepsilon ^2} + \frac{1}{\varepsilon ^2} \varrho _{\varepsilon }e(\varrho _{\varepsilon }, \vartheta _{\varepsilon }) + \frac{1}{\varepsilon ^2} \varrho _{\varepsilon }|s(\varrho _{\varepsilon }, \vartheta _{\varepsilon })| + \varrho _{\varepsilon }|\textbf{u}_\varepsilon |^2 \right) _{\textrm{res}}. \end{aligned}$$

(5.2)

5.1 Energy estimates

In agreement with hypothesis (4.2), we have

$$\begin{aligned} \int _{\Omega } { E_\varepsilon \left( \varrho _{\varepsilon }(0, \cdot ), \vartheta _{\varepsilon }(0, \cdot ), \textbf{u}_\varepsilon (0, \cdot ) \Big | \overline{\varrho }, \overline{\vartheta } + \varepsilon \Theta _B , 0 \right) } \ \,\textrm{d} {x} \lesssim 1 \ \text{ independently } \text{ of }\ \varepsilon \rightarrow 0, \end{aligned}$$

(5.3)

where $\Theta _B = \Theta _B(x)$ is a suitable extension of the temperature boundary condition inside $\Omega$. Plugging this ansatz in the relative energy inequality (2.6) we obtain

$$\begin{aligned}&\left[ \int _{\Omega } { E_\varepsilon \left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon \Big | \overline{\varrho }, \overline{\vartheta } + \varepsilon \Theta _B, 0 \right) } \ \,\textrm{d} {x} \right] _{t = 0}^{t = \tau } \nonumber \\&\qquad + \int _0^\tau \int _{\Omega } { \frac{\overline{\vartheta } + \varepsilon \Theta _B}{\vartheta _{\varepsilon }} \left( {\mathbb {S}} (\vartheta _{\varepsilon }, {\mathbb {D}}_x\textbf{u}_\varepsilon ) : {\mathbb {D}}_x\textbf{u}_\varepsilon + \frac{1}{\varepsilon ^2} \frac{\kappa (\vartheta _{\varepsilon }) |\nabla _{x}\vartheta _{\varepsilon }|^2 }{\vartheta _{\varepsilon }} \right) } \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\quad \le - \frac{1}{\varepsilon } \int _0^\tau \int _{\Omega } { \left( \varrho _{\varepsilon }(s(\varrho _{\varepsilon }, \vartheta _{\varepsilon }) - s(\overline{\varrho }, \overline{\vartheta } + \varepsilon {\Theta _B})) \textbf{u}_\varepsilon \cdot \nabla _{x}\Theta _B - \frac{\kappa (\vartheta _{\varepsilon }) \nabla _{x}\vartheta _{\varepsilon }}{\vartheta _{\varepsilon }} \cdot \nabla _{x}\Theta _B \right) } \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\qquad + \int _0^\tau \int _{\Omega } { \varrho _{\varepsilon }\frac{1}{\varepsilon } \nabla _{x}G \cdot \textbf{u}_\varepsilon } \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\qquad - \frac{1}{\varepsilon } \int _0^\tau \int _{\Omega } { \left( \frac{ \partial p(\overline{\varrho }, \overline{\vartheta } + \varepsilon \Theta _B )}{\partial \vartheta } - \frac{ \partial p(\overline{\varrho }, \overline{\vartheta } )}{\partial \vartheta } \right) \frac{\varrho _{\varepsilon }}{\overline{\varrho }} \textbf{u}_\varepsilon \cdot \nabla _{x}\Theta _B } \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\qquad - \frac{1}{\varepsilon } \int _0^\tau \int _{\Omega } { \frac{ \partial p(\overline{\varrho }, \overline{\vartheta } )}{\partial \vartheta } \frac{\varrho _{\varepsilon }}{\overline{\varrho }} \textbf{u}_\varepsilon \cdot \nabla _{x}\Theta _B } \ \,\textrm{d} {x} \,\textrm{d} t. \end{aligned}$$

(5.4)

Our goal is to control the integrals on the right-hand side to apply Grönwall’s argument. To this end, we fix the compact set K determining the essential and residual component to contain the point $(\overline{\varrho }, \overline{\vartheta })$ in its interior. In particular, the same is true for the range of the function $(\overline{\varrho }, \overline{\vartheta } + \varepsilon \Theta _B)$ as soon as $\varepsilon > 0$ is small enough. Accordingly, we will systematically use the coercivity of the relative energy $E_\varepsilon$ stated in (5.1), (5.2) in the estimates below. In particular, we have the estimate

$$\begin{aligned} \left| \left[ B(\varrho _{\varepsilon }, \vartheta _{\varepsilon }) - B(\overline{\varrho }, \overline{\vartheta } + \varepsilon \Theta _B) \right] _{\textrm{ess}} \right| ^2&\lesssim \left[ E_\varepsilon \left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon \Big | \overline{\varrho }, \overline{\vartheta } + \varepsilon \Theta _B, 0 \right) \right] _{\textrm{ess}} \\&\le E_\varepsilon \left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon \Big | \overline{\varrho }, \overline{\vartheta } + \varepsilon \Theta _B, 0 \right) \end{aligned}$$

for any $B = B(\varrho , \vartheta )$ locally Lipschitz in $(0, \infty )^2$.

5.1.1 Estimates

Step 1: First,

$$\begin{aligned}&\frac{1}{\varepsilon } \int _{\Omega } { \left| \varrho _{\varepsilon }(s(\varrho _{\varepsilon }, \vartheta _{\varepsilon }) - s(\overline{\varrho }, \overline{\vartheta } + \varepsilon {\Theta _B})) \textbf{u}_\varepsilon \cdot \nabla _{x}\Theta _B \right| } \ \,\textrm{d} {x} \\&\quad \lesssim \frac{1}{\varepsilon } \int _{\Omega } { \left| \left[ \varrho _{\varepsilon }(s(\varrho _{\varepsilon }, \vartheta _{\varepsilon }) - s(\overline{\varrho }, \overline{\vartheta } + \varepsilon {\Theta _B})) \textbf{u}_\varepsilon \right] _{\textrm{ess}} \right| } \ \,\textrm{d} {x} \\&\qquad + \frac{1}{\varepsilon } \int _{\Omega } { \left| \left[ \varrho _{\varepsilon }(s(\varrho _{\varepsilon }, \vartheta _{\varepsilon }) - s(\overline{\varrho }, \overline{\vartheta } + \varepsilon {\Theta _B})) \textbf{u}_\varepsilon \right] _{\textrm{res}} \right| } \ \,\textrm{d} {x}, \end{aligned}$$

where

$$\begin{aligned}&\frac{1}{\varepsilon } \int _{\Omega } { \left| \left[ \varrho _{\varepsilon }(s(\varrho _{\varepsilon }, \vartheta _{\varepsilon }) - s(\overline{\varrho }, \overline{\vartheta } + \varepsilon {\Theta _B})) \textbf{u}_\varepsilon \right] _{\textrm{ess}} \right| } \ \,\textrm{d} {x} \nonumber \\&\quad \lesssim \frac{1}{\varepsilon ^2} \int _{\Omega } { \left| \left[ (s(\varrho _{\varepsilon }, \vartheta _{\varepsilon }) - s(\overline{\varrho }, \overline{\vartheta } + \varepsilon {\Theta _B})) \right] _{\textrm{ess}} \right| ^2 } \ \,\textrm{d} {x} \nonumber \\&\qquad + \int _{\Omega } { \varrho _{\varepsilon }|\textbf{u}_\varepsilon |^2 } \ \,\textrm{d} {x} \lesssim \int _{\Omega } { E_\varepsilon \left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon \Big | \overline{\varrho }, \overline{\vartheta } + \varepsilon \Theta _B, 0 \right) } \ \,\textrm{d} {x}, \end{aligned}$$

(5.5)

and

$$\begin{aligned}&\frac{1}{\varepsilon } \int _{\Omega } { \left| \left[ \varrho _{\varepsilon }(s(\varrho _{\varepsilon }, \vartheta _{\varepsilon }) - s(\overline{\varrho }, \overline{\vartheta } + \varepsilon {\Theta _B})) \textbf{u}_\varepsilon \right] _{\textrm{res}} \right| } \ \,\textrm{d} {x}\\&\quad \lesssim \frac{1}{\varepsilon } \int _{\Omega } { \left[ \varrho _{\varepsilon }|\textbf{u}_\varepsilon | \right] _{\textrm{res}} } \ \,\textrm{d} {x} + \frac{1}{\varepsilon } \int _{\Omega } { \left[ \varrho _{\varepsilon }s_{\textrm{m}}(\varrho _{\varepsilon }, \vartheta _{\varepsilon }) |\textbf{u}_\varepsilon | \right] _{\textrm{res}} } \ \,\textrm{d} {x} + \frac{1}{\varepsilon } \int _{\Omega } { \left[ \vartheta _{\varepsilon }^3 |\textbf{u}_\varepsilon | \right] _{\textrm{res}} } \ \,\textrm{d} {x}. \end{aligned}$$

Furthermore,

$$\begin{aligned} \frac{1}{\varepsilon } \int _{\Omega } { \left[ \varrho _{\varepsilon }|{\textbf {u}}_\varepsilon | \right] _{\text {res}} } \ \,\text {d} {x}&\lesssim \frac{1}{\varepsilon ^2} \int _{\Omega } { [\varrho _\varepsilon ]_{\text {res}} } \ \,\text {d} {x} + \int _{\Omega } {\varrho _{\varepsilon }|{\textbf {u}}_\varepsilon |^2} \ \,\text {d} {x} \nonumber \\ {}&\lesssim \int _{\Omega } { E_\varepsilon \left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, {\textbf {u}}_\varepsilon \Big | \overline{\varrho }, \overline{\vartheta } + \varepsilon \Theta _B, 0 \right) } \ \,\text {d} {x}. \end{aligned}$$

(5.6)

In view of the bounds (2.17), (2.18),

$$\begin{aligned} \frac{1}{\varepsilon } \int _{\Omega } { \left[ \varrho _{\varepsilon }s_{\textrm{m}}(\varrho _{\varepsilon }, \vartheta _{\varepsilon }) |\textbf{u}_\varepsilon | \right] _{\textrm{res}} } \ \,\textrm{d} {x}&\lesssim \frac{1}{\varepsilon ^2} \int _{\Omega } { \left[ \varrho _{\varepsilon }s^2_{\textrm{m}}(\varrho _{\varepsilon }, \vartheta _{\varepsilon }) \right] _{\textrm{res}} } \ \,\textrm{d} {x} + \int _{\Omega } {\varrho _{\varepsilon }|\textbf{u}_\varepsilon |^2} \ \,\textrm{d} {x} \nonumber \\&\lesssim \int _{\Omega } { E_\varepsilon \left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon \Big | \overline{\varrho }, \overline{\vartheta } + \varepsilon \Theta _B, 0 \right) } \ \,\textrm{d} {x}. \end{aligned}$$

(5.7)

Finally,

$$\begin{aligned} \frac{1}{\varepsilon } \int _{\Omega } { \left[ \vartheta _{\varepsilon }^3 |\textbf{u}_\varepsilon | \right] _{\textrm{res}} } \ \,\textrm{d} {x} \lesssim \delta \Vert \textbf{u}_\varepsilon \Vert _{W^{1,2}(\Omega ; {\mathbb {R}}^d)}^2 + \frac{ C(\delta ) }{\varepsilon ^2} \int _{\Omega } { [\vartheta _{\varepsilon }^6]_{\textrm{res}} } \ \,\textrm{d} {x} \end{aligned}$$

for any $\delta > 0$. Thus if $\delta > 0$ is chosen small enough, the first integral is controlled by the viscosity dissipation on the left-hand side of (5.4). Next, in accordance with hypothesis (2.16),

$$\begin{aligned} \int _{\Omega } { \frac{\kappa (\vartheta _{\varepsilon }) |\nabla _{x}\vartheta _{\varepsilon }|^2 }{\vartheta _{\varepsilon }^2} } \ \,\text {d} {x} \gtrsim \int _{\Omega } { |\nabla _{x}\log (\vartheta _{\varepsilon }) |^2 + |\nabla _{x}\vartheta _{\varepsilon }^{\frac{\beta }{2}} |^2 } \ \,\text {d} {x}, \ \beta > 6. \end{aligned}$$

(5.8)

Consequently, as the measure of the residual set is controlled by the relative energy (cf. (5.2)), we get

$$\begin{aligned} \frac{1}{\varepsilon ^2}&\int _{\Omega } { [\vartheta _{\varepsilon }^6]_{\textrm{res}} } \ \,\textrm{d} {x} \lesssim \frac{1}{\varepsilon ^2} \int _{\Omega } { [\vartheta _{\varepsilon }^3 - \overline{\vartheta }^3 ]_{\textrm{res}}^2 } \ \,\textrm{d} {x} + \frac{1}{\varepsilon ^2} \int _{\Omega } { [\overline{\vartheta }^3]_{\textrm{res}}^2 } \ \,\textrm{d} {x} \nonumber \\&\lesssim \frac{1}{\varepsilon ^2} \int _{\Omega } { |\nabla _{x}\vartheta _{\varepsilon }^{3} |^2 } \ \,\textrm{d} {x} + \frac{1}{\varepsilon ^2} \int _{\Omega } [\vartheta _{\varepsilon }^3 - \overline{\vartheta }^3 ]_{\textrm{ess}}^2 + \frac{1}{\varepsilon ^2} \int _{\Omega } { [\overline{\vartheta }^3]_{\textrm{res}}^2 } \ \,\textrm{d} {x}\nonumber \\&\lesssim \frac{\delta }{\varepsilon ^2} \int _{\Omega } { \frac{\kappa (\vartheta _{\varepsilon }) | \nabla _{x}\vartheta _{\varepsilon }|^2 }{\vartheta _{\varepsilon }^2} } \ \,\textrm{d} {x} + C(\delta )\int _{\Omega } { E_\varepsilon \left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon \Big | \overline{\varrho }, \overline{\vartheta } + \varepsilon \Theta _B, 0 \right) } \ \,\textrm{d} {x} \end{aligned}$$

(5.9)

for any $\delta > 0$.

Step 2: In accordance with hypothesis (2.16),

$$\begin{aligned} \frac{1}{\varepsilon } \left| \int _{\Omega } { \frac{\kappa (\vartheta _{\varepsilon }) }{\vartheta _{\varepsilon }} \nabla _{x}\vartheta _{\varepsilon }} \ \,\textrm{d} {x} \right| \lesssim \frac{1}{\varepsilon } \int _{\Omega } { |\nabla _{x}(\log (\vartheta _{\varepsilon }))| + \vartheta _{\varepsilon }^{\beta - 1} |\nabla _{x}\vartheta _{\varepsilon }| } \ \,\textrm{d} {x}, \end{aligned}$$

where

$$\begin{aligned} \frac{1}{\varepsilon } \int _{\Omega } { \left| \nabla _{x}(\log (\vartheta _{\varepsilon })) \right| } \ \,\textrm{d} {x} \lesssim \frac{\delta }{\varepsilon ^2} \int _{\Omega } { \left| \nabla _{x}(\log (\vartheta _{\varepsilon })) \right| ^2 } \ \,\textrm{d} {x} + C(\delta ) \end{aligned}$$

(5.10)

for any $\delta > 0$; hence the integral is controlled by dissipation.

Next,

$$\begin{aligned} \frac{1}{\varepsilon } \int _{\Omega } { \vartheta _{\varepsilon }^{\beta - 1}{\nabla _{x}\vartheta _{\varepsilon }} } \ \,\textrm{d} {x} =\frac{1}{\varepsilon } \int _{\Omega } { \vartheta ^{\frac{\beta }{2}} \nabla _{x}\vartheta _{\varepsilon }^{\frac{\beta }{2}} } \ \,\textrm{d} {x} \le \frac{\delta }{\varepsilon ^2} \int _{\Omega } { |\nabla _{x}\vartheta _{\varepsilon }^{\frac{\beta }{2}} |^2 } \ \,\textrm{d} {x} + C(\delta ) \int _{\Omega } { |\vartheta _{\varepsilon }^{\frac{\beta }{2} }|^2 } \ \,\textrm{d} {x}, \end{aligned}$$

(5.11)

where the first term is controlled by dissipation and the second one by Poincaré’s inequality

$$\begin{aligned} \int _{\Omega } { |\vartheta _{\varepsilon }^{\frac{\beta }{2} } |^2 } \ \,\textrm{d} {x} \lesssim \int _{\Omega } { |\nabla _{x}\vartheta _{\varepsilon }^{\frac{\beta }{2}} |^2 } \ \,\textrm{d} {x} + \int _{\partial \Omega } (\overline{\vartheta } + \varepsilon \Theta _B )^\beta \ \textrm{d}\sigma _x. \end{aligned}$$

(5.12)

Step 3: We have

$$\begin{aligned} \frac{1}{\varepsilon } \int _0^\tau \int _{\Omega } { \varrho _{\varepsilon }\nabla _{x}G \cdot \textbf{u}_\varepsilon } \ \,\textrm{d} {x} \,\textrm{d} t&= - \frac{1}{\varepsilon } \int _0^\tau \int _{\Omega } { G \textrm{div}_{x}(\varrho _{\varepsilon }\textbf{u}_\varepsilon ) } \ \,\textrm{d} {x} \,\textrm{d} t\\&= {\int _0^\tau \partial _t \int _{\Omega } { \frac{1}{\varepsilon } (\varrho _{\varepsilon }- \overline{\varrho }) G } \ \,\textrm{d} {x} \,\textrm{d} t}\\&= \left[ \int _{\Omega } { \frac{\varrho _{\varepsilon }- \overline{\varrho }}{\varepsilon } G } \ \,\textrm{d} {x} \right] _{t = 0}^{t=\tau }. \end{aligned}$$

Seeing that

$$\begin{aligned} E_\varepsilon \left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon \Big | \overline{\varrho }, \overline{\vartheta } +\varepsilon \Theta _B, 0 \right) - c_1&\lesssim E_\varepsilon \left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon \Big | \overline{\varrho }, \overline{\vartheta } + \varepsilon \Theta _B, 0 \right) - \int _{\Omega } { \frac{\varrho _{\varepsilon }- \overline{\varrho }}{\varepsilon } G } \ \,\textrm{d} {x} \nonumber \\&\lesssim E_\varepsilon \left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon \Big | \overline{\varrho }, \overline{\vartheta } + \varepsilon \Theta _B, 0 \right) + 1 \end{aligned}$$

(5.13)

we can add this term to the relative energy on the left-hand side of (2.6).

Step 4:

$$\begin{aligned}&\frac{1}{\varepsilon } \left| \int _{\Omega } { \left( \frac{ \partial p(\overline{\varrho }, \overline{\vartheta } + \varepsilon \Theta _B )}{\partial \vartheta } - \frac{ \partial p(\overline{\varrho }, \overline{\vartheta } )}{\partial \vartheta } \right) \frac{\varrho _{\varepsilon }}{\overline{\varrho }} \textbf{u}_\varepsilon \cdot \nabla _{x}\Theta _B } \ \,\textrm{d} {x} \right| \nonumber \\&\lesssim \int _{\Omega } { \varrho _{\varepsilon }|\textbf{u}_\varepsilon | } \ \,\textrm{d} {x} \lesssim \int _{\Omega } { \varrho _{\varepsilon }} \ \,\textrm{d} {x} + \int _{\Omega } { \varrho _{\varepsilon }|\textbf{u}_\varepsilon |^2 } \ \,\textrm{d} {x} \lesssim \int _{\Omega } { E_\varepsilon \left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon \Big | \overline{\varrho }, \overline{\vartheta } + \varepsilon \Theta _B, 0 \right) } \ \,\textrm{d} {x} + 1. \end{aligned}$$

(5.14)

Step 5: The last integral on the right-hand side of (5.4) can be handled exactly as in Step 3.

5.1.2 Conclusion, uniform bounds

In view of the estimates obtained in the previous section, we may apply Grönwall’s lemma to the relative energy inequality (5.4). As the initial data satisfy (5.3), we deduce the following bounds independent of the scaling parameter $\varepsilon \rightarrow 0$:

$$\begin{aligned} \textrm{ess} \sup _{t \in (0,T)} \int _{\Omega } { E_\varepsilon \left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon \Big | \overline{\varrho }, \overline{\vartheta } + \varepsilon \Theta _B, 0 \right) } \ \,\textrm{d} {x}&\lesssim 1, \end{aligned}$$

(5.15)

$$\begin{aligned} \int _0^T \Vert \textbf{u}_\varepsilon \Vert ^2_{W^{1,2}_0 (\Omega ; {\mathbb {R}}^d) } \,\textrm{d} t&\lesssim 1, \end{aligned}$$

(5.16)

$$\begin{aligned} \frac{1}{\varepsilon ^2} \int _0^T \left( \Vert \nabla _{x}\log (\vartheta _{\varepsilon }) \Vert ^2_{L^2(\Omega ; {\mathbb {R}}^d)} + \Vert \nabla _{x}\vartheta _{\varepsilon }^{\frac{\beta }{2}} \Vert ^2_{L^2(\Omega ; {\mathbb {R}}^d)} \right)&\lesssim 1. \end{aligned}$$

(5.17)

Next, it follows from (5.15) that the measure of the residual set shrinks to zero, specifically

$$\begin{aligned} \frac{1}{\varepsilon ^2} \textrm{ess} \sup _{t \in (0,T)} \int _{\Omega } { [1]_{\textrm{res}}} \ \,\textrm{d} {x} \lesssim 1. \end{aligned}$$

(5.18)

In addition, we get from (5.15):

$$\begin{aligned} \textrm{ess} \sup _{t \in (0,T)} \int _{\Omega } { \varrho _{\varepsilon }|\textbf{u}_\varepsilon |^2 } \ \,\textrm{d} {x}&\lesssim 1, \nonumber \\ \textrm{ess} \sup _{t \in (0,T)} \left\| \left[ \frac{\varrho _{\varepsilon }- \overline{\varrho }}{\varepsilon } \right] _{\textrm{ess}} \right\| _{L^2(\Omega )}&\lesssim 1,\nonumber \\ \textrm{ess} \sup _{t \in (0,T)} \left\| \left[ \frac{\vartheta _{\varepsilon }- \overline{\vartheta }}{\varepsilon } \right] _{\textrm{ess}} \right\| _{L^2(\Omega )}&\lesssim 1,\nonumber \\ \frac{1}{\varepsilon ^2} \textrm{ess} \sup _{t \in (0,T)} \Vert [\varrho _{\varepsilon }]_{\textrm{res}} \Vert ^{\frac{5}{3}}_{L^{\frac{5}{3}}(\Omega )} + \frac{1}{\varepsilon ^2} \textrm{ess} \sup _{t \in (0,T)} \Vert [\vartheta _{\varepsilon }]_{\textrm{res}} \Vert ^{4}_{L^{4}(\Omega )}&\lesssim 1. \end{aligned}$$

(5.19)

Combining (5.17), (5.18), and (5.19), we conclude

$$\begin{aligned} \int _0^T \left\| \frac{\log (\vartheta _{\varepsilon }) - \log (\overline{\vartheta })}{\varepsilon } \right\| ^2_{W^{1,2}(\Omega )} \,\textrm{d} t+ \int _0^T \left\| \frac{ \vartheta _{\varepsilon }-\overline{\vartheta } }{\varepsilon } \right\| ^2_{W^{1,2}(\Omega )} \,\textrm{d} t\lesssim 1. \end{aligned}$$

(5.20)

Finally, we claim the bound on the entropy flux

$$\begin{aligned} \int _0^T \left\| \left[ \frac{\kappa (\vartheta _{\varepsilon }) }{\vartheta _{\varepsilon }} \right] _{\textrm{res}} \frac{\nabla _{x}\vartheta _{\varepsilon }}{\varepsilon } \right\| ^q_{L^q(\Omega ; {\mathbb {R}}^d)} \,\textrm{d} t\lesssim 1 \ \text{ for } \text{ some }\ q > 1. \end{aligned}$$

(5.21)

Indeed we have

$$\begin{aligned} \left| \left[ \frac{\kappa (\vartheta _{\varepsilon }) }{\vartheta _{\varepsilon }} \right] _{\textrm{res}} \frac{\nabla _{x}\vartheta _{\varepsilon }}{\varepsilon } \right| \lesssim \frac{1}{\varepsilon } \left| \nabla _{x}\log (\vartheta _{\varepsilon }) \right| + \frac{1}{\varepsilon } \left| \left[ \vartheta _{\varepsilon }^{\frac{\beta }{2}} \nabla _{x}\vartheta _{\varepsilon }^{\frac{\beta }{2}} \right] _{\textrm{res}} \right| , \end{aligned}$$

where the former term on the right-hand side is controlled via (5.20). As for the latter, we deduce from (5.17) that

$$\begin{aligned} \left\| \frac{1}{\varepsilon } \nabla _{x}\vartheta _{\varepsilon }^{\frac{\beta }{2}} \right\| _{L^2((0,T) \times \Omega ; {\mathbb {R}}^d)} \lesssim 1; \end{aligned}$$

hence it is enough to check

$$\begin{aligned} \left\| \left[ \vartheta _{\varepsilon }^{\frac{\beta }{2}} \right] _{\textrm{res}} \right\| _{L^r ((0,T) \times \Omega )} \lesssim 1 \ \text{ for } \text{ some }\ r > 2. \end{aligned}$$

(5.22)

To see (5.22), first observe that

$$\begin{aligned} \textrm{ess} \sup _{t \in (0,T)} \Vert [\vartheta _{\varepsilon }]_{\textrm{res}} \Vert _{L^4(\Omega )} \lesssim 1, \end{aligned}$$

and, in view of (5.17) and Poincaré inequality,

$$\begin{aligned} \left\| \vartheta _{\varepsilon }^{\frac{\beta }{2}} \right\| _{L^2(0,T; L^6(\Omega ))} \lesssim 1 \ (\text{ for }\ d = 3). \end{aligned}$$

Consequently, (5.22) follows by interpolation.

6 Convergence to the Target System

Our ultimate goal is to perform the limit $\varepsilon \rightarrow 0$. We proceed in two steps.

6.1 Weak convergence

In view of the uniform bounds established in Sect. 5.1.2, we may infer

$$\begin{aligned} \varrho _{\varepsilon }&\rightarrow \overline{\varrho } \ \text{ in }\ L^{\frac{5}{3}}(\Omega ) \ \text{ uniformly } \text{ for }\ t \in (0,T), \end{aligned}$$

(6.1)

$$\begin{aligned} \vartheta _{\varepsilon }&\rightarrow \overline{\vartheta } \ \text{ in }\ L^2(0,T; W^{1,2}(\Omega )), \end{aligned}$$

(6.2)

$$\begin{aligned} \textbf{u}_\varepsilon&\rightarrow \textbf{u}\ \text{ weakly } \text{ in }\ L^2(0,T; W^{1,2}_0(\Omega ; {\mathbb {R}}^d)), \end{aligned}$$

(6.3)

where (6.3) may require extraction of a suitable subsequence. As we shall eventually see, the limit velocity $\textbf{u}= \textbf{U}$ is unique so that the convergence is, in fact, unconditional. In addition, we may let $\varepsilon \rightarrow 0$ in the weak formulation of the equation of continuity (2.2) to deduce

$$\begin{aligned} \textrm{div}_{x}\textbf{u}= 0. \end{aligned}$$

(6.4)

Next, we use (5.19), (5.20) to obtain (a priori for suitable subsequences),

$$\begin{aligned} \frac{ \varrho _{\varepsilon }- \overline{\varrho } }{\varepsilon }&= \left[ \frac{ \varrho _{\varepsilon }- \overline{\varrho } }{\varepsilon } \right] _{\textrm{ess}} +\left[ \frac{ \varrho _{\varepsilon }- \overline{\varrho } }{\varepsilon } \right] _{\textrm{res}}, \nonumber \\ \left[ \frac{ \varrho _{\varepsilon }- \overline{\varrho } }{\varepsilon } \right] _{\textrm{ess}}&\rightarrow \mathfrak {R} \ \text{ weakly-(*) } \text{ in } \ L^\infty (0,T; L^2(\Omega )),\nonumber \\ \left[ \frac{ \varrho _{\varepsilon }- \overline{\varrho } }{\varepsilon } \right] _{\textrm{res}}&\rightarrow 0 \ \text{ in }\ L^\infty (0,T; L^{\frac{5}{3}}(\Omega )), \end{aligned}$$

(6.5)

$$\begin{aligned} \frac{ \vartheta _{\varepsilon }- \overline{\vartheta } }{\varepsilon }&\rightarrow \mathfrak {T} \ \text{ weakly } \text{ in }\ L^2(0,T; W^{1,2}(\Omega )) \ \text{ and } \text{ weakly-(*) } \text{ in }\ L^\infty (0,T; L^2(\Omega )). \end{aligned}$$

(6.6)

Moreover, in view of (4.1),

$$\begin{aligned} \mathfrak {T}|_{\partial \Omega } = \Theta _B. \end{aligned}$$

(6.7)

Finally, we perform the limit in the rescaled momentum equation (1.2) to deduce

$$\begin{aligned} \frac{\partial p(\overline{\varrho }, \overline{\vartheta }) }{\partial \varrho } \nabla _{x}\mathfrak {R} +\frac{\partial p(\overline{\varrho }, \overline{\vartheta }) }{\partial \vartheta } \nabla _{x}\mathfrak {T} = \overline{\varrho } \nabla _{x}G \end{aligned}$$

(6.8)

in the sense of distributions. In particular, it follows from (6.8) that

$$\begin{aligned} \mathfrak {R} \in L^2(0,T; W^{1,2}(\Omega )). \end{aligned}$$

(6.9)

6.2 Strong convergence

First, it is more convenient to rewrite the target OB system in terms of the variable

Accordingly, we get

(6.10)

together with the Boussinesq relation

$$\begin{aligned} \frac{\partial p(\overline{\varrho }, \overline{\vartheta } ) }{\partial \varrho } \nabla _{x}r +\frac{\partial p(\overline{\varrho }, \overline{\vartheta } ) }{\partial \vartheta } \nabla _{x}{\mathcal {T}}= \overline{\varrho } \nabla _{x}G,\ \int _{\Omega } { r } \ \,\textrm{d} {x} = 0, \end{aligned}$$

(6.11)

the boundary conditions

$$\begin{aligned} \textbf{U}|_{\partial \Omega } = 0,\ {\mathcal {T}}|_{\partial \Omega } = \Theta _B, \end{aligned}$$

(6.12)

and the initial conditions

$$\begin{aligned} \textbf{U}(0, \cdot ) = \textbf{U}_0,\ {\mathcal {T}}(0, \cdot ) = \mathfrak {T}_0. \end{aligned}$$

(6.13)

In accordance with Proposition 3.1 and hypotheses (4.4), (4.5), the problem (6.10)–(6.13) admits a unique regular solution on a time interval $[0, T_{\textrm{max}})$, where $T_{\textrm{max}} > 0$ and $T_{\textrm{max}} = \infty$ if $d=2$.

6.3 Relative energy

To complete the proof of Theorem 4.1, we use the relative energy inequality (2.6), with the ansatz

$$\begin{aligned} E_\varepsilon \left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon \ \Big | \overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}, \textbf{U}\right) . \end{aligned}$$

In accordance with our choice of the initial data,

$$\begin{aligned} \int _{\Omega } { E_\varepsilon \left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon \ \Big | \overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}, \textbf{U}\right) (0, \cdot ) } \ \,\textrm{d} {x} \rightarrow 0 \ \text{ as }\ \varepsilon \rightarrow 0. \end{aligned}$$

(6.14)

Our goal is to show that $\displaystyle \int _{\Omega } {E_\varepsilon \left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon \ \Big | \overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}, \textbf{U}\right) (\tau ,\cdot )} \ \,\textrm{d} {x}\rightarrow 0$, which finally yields

$$\begin{aligned} {\mathcal {T}}=\lim _{\varepsilon \rightarrow 0} \frac{\vartheta _{\varepsilon }-\overline{\vartheta }}{\varepsilon } =\mathfrak {T},\ r=\lim _{\varepsilon \rightarrow 0} \frac{\varrho _{\varepsilon }-\overline{\varrho }}{\varepsilon } =\mathfrak {R},\ \lim _{\varepsilon \rightarrow 0} \textbf{u}_\varepsilon =\textbf{U}, \end{aligned}$$

and the trio $(\mathfrak {R},\mathfrak {T},\textbf{U})=(r, {\mathcal {T}}, \textbf{U})$ is the strong solution to the OB system (1.25)–(1.27).

Step 1: Plugging our ansatz in the relative energy inequality (2.6) and using $\textrm{div}_{x}\textbf{U}= 0$ we get

$$\begin{aligned}&\left[ \int _{\Omega } { E_\varepsilon \left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon \Big | \overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}, \textbf{U}\right) } \ \,\textrm{d} {x} \right] _{t = 0}^{t = \tau }\nonumber \\&\qquad + \int _0^\tau \int _{\Omega } { \frac{\overline{\vartheta } + \varepsilon {\mathcal {T}}}{\vartheta _{\varepsilon }} \left( {\mathbb {S}} (\vartheta _{\varepsilon }, {\mathbb {D}}_x\textbf{u}_\varepsilon ) : {\mathbb {D}}_x\textbf{u}_\varepsilon + \frac{1}{\varepsilon ^2} \frac{\kappa (\vartheta _{\varepsilon }) \nabla _{x}\vartheta _{\varepsilon }\cdot \nabla _{x}\vartheta _{\varepsilon }}{\vartheta _{\varepsilon }} \right) } \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\quad \le - \frac{1}{\varepsilon } \int _0^\tau \int _{\Omega } { \varrho _{\varepsilon }\Big [ (s(\varrho _{\varepsilon }, \vartheta _{\varepsilon }) - s( \overline{\varrho } + \varepsilon r,\overline{\vartheta } + \varepsilon {\mathcal {T}}) \Big ] \partial _t {\mathcal {T}}} \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\qquad - \frac{1}{\varepsilon } \int _0^\tau \int _{\Omega } { \varrho _{\varepsilon }\Big [ s(\varrho _{\varepsilon }, \vartheta _{\varepsilon }) - s( \overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}) \Big ] \textbf{u}_\varepsilon \cdot \nabla _{x}{\mathcal {T}}} \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\qquad + \frac{1}{\varepsilon } \int _0^\tau \int _{\Omega } {\frac{\kappa (\vartheta _{\varepsilon })}{\vartheta _{\varepsilon }} \nabla _{x}\vartheta _{\varepsilon }\cdot \nabla _{x}{\mathcal {T}}} \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\qquad - \int _0^\tau \int _{\Omega } { \Big [ \varrho _{\varepsilon }(\textbf{u}_\varepsilon - \textbf{U}) \otimes (\textbf{u}_\varepsilon - \textbf{U}) - {\mathbb {S}}(\vartheta _{\varepsilon }, {\mathbb {D}}_x\textbf{u}_\varepsilon ) \Big ] : {\mathbb {D}}_x\textbf{U}} \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\qquad + \int _0^\tau \int _{\Omega } { \varrho _{\varepsilon }\left[ \frac{1}{\varepsilon } \nabla _{x}G - \partial _t \textbf{U}- (\textbf{U}\cdot \nabla _{x}) \textbf{U}\right] \cdot (\textbf{u}_\varepsilon - \textbf{U}) } \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\qquad + \frac{1}{\varepsilon ^2} \int _0^\tau \int _{\Omega } \left[ \left( 1 - \frac{\varrho _{\varepsilon }}{\overline{\varrho } + \varepsilon r} \right) \partial _t p ( \overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}) \right. \nonumber \\&\qquad \qquad \qquad \qquad \qquad \left. - \frac{\varrho _{\varepsilon }}{\overline{\varrho } + \varepsilon r} \textbf{u}_\varepsilon \cdot \nabla _{x}p( \overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}) \right] \ \,\textrm{d} {x}\,\textrm{d} t. \end{aligned}$$

(6.15)

Step 2: As r, $\textbf{U}$ satisfy the momentum equation

$$\begin{aligned} - \overline{\varrho } \left( \partial _t \textbf{U}+ \textbf{U}\cdot \nabla _{x}\textbf{U}\right) = \nabla _{x}\Pi - \textrm{div}_{x}{\mathbb {S}}(\overline{\vartheta }, \nabla _{x}\textbf{U}) - r \nabla _{x}G, \end{aligned}$$

we get

$$\begin{aligned}&\int _{\Omega } {\varrho _{\varepsilon }\left[ \frac{1}{\varepsilon } \nabla _{x}G - \partial _t \textbf{U}- (\textbf{U}\cdot \nabla _{x}) \textbf{U}\right] \cdot (\textbf{u}_\varepsilon - \textbf{U}) } \ \,\textrm{d} {x} \\&\quad = \int _{\Omega } { \frac{\varrho _{\varepsilon }}{\overline{\varrho }} \left[ \frac{1}{\varepsilon } \overline{\varrho } \nabla _{x}G + \nabla _{x}\Pi - \textrm{div}_{x}{\mathbb {S}}(\overline{\vartheta }, \nabla _{x}\textbf{U}) - r \nabla _{x}G \right] \cdot (\textbf{u}_\varepsilon - \textbf{U}) } \ \,\textrm{d} {x}. \end{aligned}$$

Thus we can use the convergence established in (6.1)–(6.3) to rewrite (6.15) in the form

$$\begin{aligned}&\int _{\Omega } { E_\varepsilon \left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon \Big | \overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}, \textbf{U}\right) (\tau , \cdot ) } \ \,\textrm{d} {x} \\&\qquad + \int _0^\tau \int _{\Omega } { \Big ( {\mathbb {S}} (\overline{\vartheta }, {\mathbb {D}}_x\textbf{u}_\varepsilon ) - {\mathbb {S}} (\overline{\vartheta }, {\mathbb {D}}_x\textbf{U}) \Big ) : \Big ( {\mathbb {D}}_x\textbf{u}_\varepsilon - {\mathbb {D}}_x\textbf{U}\Big ) } \ \,\textrm{d} {x} \,\textrm{d} t\\&\qquad +\int _0^\tau \int _{\Omega } { \left( \frac{\overline{\vartheta } {+} \varepsilon {\mathcal {T}}}{\vartheta _{\varepsilon }^2} \right) \frac{\kappa (\vartheta _{\varepsilon }) \nabla _{x}\vartheta _{\varepsilon }\cdot \nabla _{x}\vartheta _{\varepsilon }}{\varepsilon ^2} } \ \,\textrm{d} {x} \,\textrm{d} t{-} \int _0^\tau \int _{\Omega } { \frac{\kappa (\vartheta _{\varepsilon }) }{\vartheta _{\varepsilon }} \frac{\nabla _{x}\vartheta _{\varepsilon }}{\varepsilon } \cdot \nabla _{x}{\mathcal {T}}} \ \,\textrm{d} {x} \,\textrm{d} t\\&\quad \le - \frac{1}{\varepsilon } \int _0^\tau \int _{\Omega } { \varrho _{\varepsilon }\Big [ (s(\varrho _{\varepsilon }, \vartheta _{\varepsilon }) - s( \overline{\varrho } + \varepsilon r,\overline{\vartheta } + \varepsilon {\mathcal {T}}) \Big ] \partial _t {\mathcal {T}}} \ \,\textrm{d} {x} \,\textrm{d} t\\&\qquad - \frac{1}{\varepsilon } \int _0^\tau \int _{\Omega } { \varrho _{\varepsilon }\Big [ s(\varrho _{\varepsilon }, \vartheta _{\varepsilon }) - s( \overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}) \Big ] \textbf{u}_\varepsilon \cdot \nabla _{x}{\mathcal {T}}} \ \,\textrm{d} {x} \,\textrm{d} t\\&\qquad + \int _0^\tau \int _{\Omega } { \frac{\varrho _{\varepsilon }}{\overline{\varrho }} \nabla _{x}\Pi \cdot (\textbf{u}_\varepsilon - \textbf{U}) } \ \,\textrm{d} {x} \,\textrm{d} t\\&\qquad + \int _0^\tau \int _{\Omega } { \frac{\varrho _{\varepsilon }}{\overline{\varrho }} \left[ \frac{1}{\varepsilon } \overline{\varrho } \nabla _{x}G - r \nabla _{x}G \right] \cdot (\textbf{u}_\varepsilon - \textbf{U}) } \ \,\textrm{d} {x} \,\textrm{d} t\\&\qquad + \frac{1}{\varepsilon ^2} \int _0^\tau \int _{\Omega } \left[ \left( 1 - \frac{\varrho _{\varepsilon }}{\overline{\varrho } + \varepsilon r} \right) \partial _t p(\overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}})\right. \\&\qquad \qquad \qquad \qquad \qquad \left. - \frac{\varrho _{\varepsilon }}{\overline{\varrho } + \varepsilon r} \textbf{u}_\varepsilon \cdot \nabla _{x}p(\overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}) \right] \ \,\textrm{d} {x}\,\textrm{d} t\\&\qquad + C \int _0^\tau \int _{\Omega } {E_\varepsilon \left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon \Big | \overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}, \textbf{U}\right) } \ \,\textrm{d} {x} \,\textrm{d} t+ {\mathcal {O}}(\varepsilon ), \end{aligned}$$

where the symbol ${\mathcal {O}}(\varepsilon )$ denotes a generic error, ${\mathcal {O}}(\varepsilon ) \rightarrow 0$ as $\varepsilon \rightarrow 0$. Note that the convective term

$$\begin{aligned} \left| \int _{\Omega } { \varrho _{\varepsilon }(\textbf{u}_\varepsilon - \textbf{U}) \otimes (\textbf{u}_\varepsilon - \textbf{U}) } \ \,\textrm{d} {x} \right| \le \int _{\Omega } { \varrho _{\varepsilon }|\textbf{u}_\varepsilon - \textbf{U} |^2 } \ \,\textrm{d} {x} \end{aligned}$$

is controlled by the relative energy.

In addition, in view of the convergences (6.1)–(6.3), we conclude

$$\begin{aligned} \int _0^\tau \int _{\Omega } { \frac{\varrho _{\varepsilon }}{\overline{\varrho }} \nabla _{x}\Pi \cdot (\textbf{u}_\varepsilon - \textbf{U}) } \ \,\textrm{d} {x} \,\textrm{d} t= {\mathcal {O}}(\varepsilon ); \end{aligned}$$

hence

$$\begin{aligned}&\int _{\Omega } { E_\varepsilon \left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon \Big | \overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}, \textbf{U}\right) (\tau , \cdot ) } \ \,\textrm{d} {x} \nonumber \\&\qquad + \int _0^\tau \int _{\Omega } { \Big ( {\mathbb {S}} (\overline{\vartheta }, {\mathbb {D}}_x\textbf{u}_\varepsilon ) - {\mathbb {S}} (\overline{\vartheta }, {\mathbb {D}}_x\textbf{U}) \Big ) : \Big ( {\mathbb {D}}_x\textbf{u}_\varepsilon - {\mathbb {D}}_x\textbf{U}\Big ) } \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\qquad +\int _0^\tau \int _{\Omega } { \left( \frac{\overline{\vartheta } + \varepsilon {\mathcal {T}}}{\vartheta _{\varepsilon }^2} \right) \frac{\kappa (\vartheta _{\varepsilon }) \nabla _{x}\vartheta _{\varepsilon }\cdot \nabla _{x}\vartheta _{\varepsilon }}{\varepsilon ^2} } \ \,\textrm{d} {x} \,\textrm{d} t- \int _0^\tau \int _{\Omega } { \frac{\kappa (\vartheta _{\varepsilon }) }{\vartheta _{\varepsilon }} \frac{\nabla _{x}\vartheta _{\varepsilon }}{\varepsilon } \cdot \nabla _{x}{\mathcal {T}}} \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\quad \le - \frac{1}{\varepsilon } \int _0^\tau \int _{\Omega } { \varrho _{\varepsilon }\Big [ (s(\varrho _{\varepsilon }, \vartheta _{\varepsilon }) - s( \overline{\varrho } + \varepsilon r,\overline{\vartheta } + \varepsilon {\mathcal {T}}) \Big ] \partial _t {\mathcal {T}}} \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\qquad - \frac{1}{\varepsilon } \int _0^\tau \int _{\Omega } { \varrho _{\varepsilon }\Big [ s(\varrho _{\varepsilon }, \vartheta _{\varepsilon }) - s( \overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}) \Big ] \textbf{u}_\varepsilon \cdot \nabla _{x}{\mathcal {T}}} \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\qquad + \int _0^\tau \int _{\Omega } { \frac{\varrho _{\varepsilon }}{\overline{\varrho }} \left[ \frac{1}{\varepsilon } \overline{\varrho } \nabla _{x}G - r \nabla _{x}G \right] \cdot (\textbf{u}_\varepsilon - \textbf{U}) } \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\qquad + \frac{1}{\varepsilon ^2} \int _0^\tau \int _{\Omega } \left[ \left( 1 - \frac{\varrho _{\varepsilon }}{\overline{\varrho } + \varepsilon r} \right) \partial _t p(\overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}) \right. \nonumber \\&\qquad \qquad \qquad \qquad \qquad \left. - \frac{\varrho _{\varepsilon }}{\overline{\varrho } + \varepsilon r} \textbf{u}_\varepsilon \cdot \nabla _{x}p(\overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}) \right] \ \,\textrm{d} {x}\,\textrm{d} t\nonumber \\&\qquad + C \int _0^\tau \int _{\Omega } {E_\varepsilon \left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon \Big | \overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}, \textbf{U}\right) } \ \,\textrm{d} {x} \,\textrm{d} t+ {\mathcal {O}}(\varepsilon ). \end{aligned}$$

(6.16)

Step 3: At this stage, we use the Boussinesq relation (6.11) to obtain

$$\begin{aligned} \frac{\partial p(\overline{\varrho }, \overline{\vartheta } ) }{\partial \varrho } r + \frac{\partial p(\overline{\varrho }, \overline{\vartheta } ) }{\partial \vartheta } {\mathcal {T}}= \overline{\varrho } G + \chi (t), \end{aligned}$$

where

since $\int _{\Omega } { G } \ \,\textrm{d} {x}=0$. Consequently,

$$\begin{aligned}&\frac{1}{\varepsilon ^2} \int _{\Omega } { \left( 1 - \frac{\varrho _{\varepsilon }}{\overline{\varrho } + \varepsilon r} \right) \partial _t p (\overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}) } \ \,\textrm{d} {x} \nonumber \\&\quad =\frac{1}{\varepsilon } \int _{\Omega } { \left( 1 - \frac{\varrho _{\varepsilon }}{\overline{\varrho } + \varepsilon r} \right) \left( \frac{\partial p (\overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}) }{\partial \varrho } \partial _t r + \frac{\partial p (\overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}) }{\partial \vartheta } \partial _t {\mathcal {T}}\right) } \ \,\textrm{d} {x} \nonumber \\&\quad = \int _{\Omega } { \frac{1}{\varepsilon } \left( 1 - \frac{\varrho _{\varepsilon }}{\overline{\varrho } + \varepsilon r} \right) \left( \frac{\partial p(\overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}})}{\partial \varrho } -\frac{\partial p (\overline{\varrho } , \overline{\vartheta } ) }{\partial \varrho } \right) \partial _t r} \ \,\textrm{d} {x} \nonumber \\&\qquad + \int _{\Omega } { \frac{1}{\varepsilon } \left( 1 - \frac{\varrho _{\varepsilon }}{\overline{\varrho } + \varepsilon r} \right) \left( \frac{\partial p (\overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}) }{\partial \vartheta } -\frac{\partial p (\overline{\varrho } , \overline{\vartheta } ) }{\partial \vartheta } \right) \partial _t {\mathcal {T}}} \ \,\textrm{d} {x} \nonumber \\&\qquad + { \frac{1}{\varepsilon } \int _{\Omega } { \left( \frac{\overline{\varrho } + \varepsilon r - \varrho _{\varepsilon }}{\overline{\varrho } + \varepsilon r} \right) \partial _t \chi } \ \,\textrm{d} {x}}, \end{aligned}$$

(6.17)

where

$$\begin{aligned}&\int _{\Omega } { \frac{1}{\varepsilon } \left( 1 - \frac{\varrho _{\varepsilon }}{\overline{\varrho } + \varepsilon r} \right) \left( \frac{\partial p(\overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}})}{\partial \varrho } -\frac{\partial p (\overline{\varrho } , \overline{\vartheta } ) }{\partial \varrho } \right) \partial _t r} \ \,\textrm{d} {x}\\&\quad + \int _{\Omega } { \frac{1}{\varepsilon } \left( 1 - \frac{\varrho _{\varepsilon }}{\overline{\varrho } + \varepsilon r} \right) \left( \frac{\partial p (\overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}) }{\partial \vartheta } -\frac{\partial p (\overline{\varrho } , \overline{\vartheta } )}{\partial \vartheta } \right) \partial _t {\mathcal {T}}} \ \,\textrm{d} {x} = {\mathcal {O}}(\varepsilon ). \end{aligned}$$

Moreover,

$$\begin{aligned} {\frac{1}{\varepsilon } \frac{\overline{\varrho } + \varepsilon r - \varrho _{\varepsilon }}{\overline{\varrho } + \varepsilon r} = - \frac{\varrho _{\varepsilon }- \overline{\varrho }}{\varepsilon ( \overline{\varrho } + \varepsilon r )} + \frac{r}{\overline{\varrho } + \varepsilon r} \rightarrow \frac{1}{\overline{\varrho }} (r - \mathfrak {R} ).} \end{aligned}$$

Seeing that

$$\begin{aligned} \int _{\Omega } { r } \ \,\textrm{d} {x} = \int _{\Omega } { \mathfrak {R} } \ \,\textrm{d} {x} = 0, \end{aligned}$$

we may infer

$$\begin{aligned} \frac{1}{\varepsilon } \int _{\Omega } { \left( \frac{ \overline{\varrho } + \varepsilon r - \varrho _{\varepsilon }}{\overline{\varrho } + \varepsilon r } \right) \partial _t \chi } \ \,\textrm{d} {x} = {\mathcal {O}}(\varepsilon ). \end{aligned}$$

Similarly,

$$\begin{aligned}&- \frac{1}{\varepsilon ^2} \int _0^\tau \int _{\Omega } { \frac{\varrho _{\varepsilon }}{\overline{\varrho } + \varepsilon r} \textbf{u}_\varepsilon \cdot \nabla _{x}p(\overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}) } \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\quad = - \frac{1}{\varepsilon } \int _0^\tau \int _{\Omega } { \frac{\varrho _{\varepsilon }}{\overline{\varrho } +\varepsilon r} \textbf{u}_\varepsilon \cdot \left( \frac{\partial p (\overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}) }{\partial \varrho } \nabla _{x}r + \frac{\partial p (\overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}) }{\partial \vartheta } \nabla _{x}{\mathcal {T}}\right) } \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\quad =- \int _0^\tau \int _{\Omega } { \frac{1}{\varepsilon } \frac{\varrho _{\varepsilon }}{\overline{\varrho } + \varepsilon r} \textbf{u}_\varepsilon \cdot \nabla _{x}r \left( \frac{\partial p (\overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}) }{\partial \varrho } -\frac{\partial p (\overline{\varrho } , \overline{\vartheta } ) }{\partial \varrho } \right) } \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\qquad - \int _0^\tau \int _{\Omega } { \frac{1}{\varepsilon } \frac{\varrho _{\varepsilon }}{\overline{\varrho } + \varepsilon r} \textbf{u}_\varepsilon \cdot \nabla _{x}{\mathcal {T}}\left( \frac{\partial p (\overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}) }{\partial \vartheta } -\frac{\partial p (\overline{\varrho } , \overline{\vartheta } ) }{\partial \vartheta } \right) } \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\qquad - \frac{1}{\varepsilon } \int _0^\tau \int _{\Omega } {\frac{\varrho _{\varepsilon }}{\overline{\varrho } + \varepsilon r} \overline{\varrho } \textbf{u}_\varepsilon \cdot \nabla _{x}G } \ \,\textrm{d} {x} \,\textrm{d} t. \end{aligned}$$

(6.18)

Using (6.1), (6.3), we perform the limit in the first integral obtaining

$$\begin{aligned}&\int _0^\tau \int _{\Omega } { \frac{1}{\varepsilon } \frac{\varrho _{\varepsilon }}{\overline{\varrho } + \varepsilon r} \textbf{u}_\varepsilon \cdot \nabla _{x}r \left( \frac{\partial p(\overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}) }{\partial \varrho } -\frac{\partial p (\overline{\varrho } , \overline{\vartheta } )}{\partial \varrho } \right) } \ \,\textrm{d} {x} \,\textrm{d} t\\&\qquad + \int _0^\tau \int _{\Omega } { \frac{1}{\varepsilon } \frac{\varrho _{\varepsilon }}{\overline{\varrho } + \varepsilon r} \textbf{u}_\varepsilon \cdot \nabla _{x}{\mathcal {T}}\left( \frac{\partial p (\overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}) }{\partial \vartheta } -\frac{\partial p (\overline{\varrho }, \overline{\vartheta } ) }{\partial \vartheta } \right) } \ \,\textrm{d} {x} \,\textrm{d} t\\&\quad = \int _0^\tau \int _{\Omega } { \textbf{u}\cdot \left( \frac{\partial ^2 p (\overline{\varrho } , \overline{\vartheta }) }{\partial ^2 \varrho } r \nabla _{x}r + \frac{\partial ^2 p (\overline{\varrho } , \overline{\vartheta }) }{\partial \varrho \partial \vartheta } \nabla _{x}(r {\mathcal {T}}) +\frac{\partial ^2 p (\overline{\varrho } , \overline{\vartheta }) }{\partial ^2 \vartheta } {{\mathcal {T}}} \nabla _{x}{\mathcal {T}}\right) } \ \,\textrm{d} {x} + {\mathcal {O}}(\varepsilon ) \\&\quad = {\mathcal {O}}(\varepsilon ) \end{aligned}$$

as $\textrm{div}_{x}\textbf{u}= 0$.

Finally,

$$\begin{aligned}&- \frac{1}{\varepsilon } \int _0^\tau \int _{\Omega } {\frac{\varrho _{\varepsilon }}{\overline{\varrho } + \varepsilon r} \overline{\varrho } \textbf{u}_\varepsilon \cdot \nabla _{x}G } \ \,\textrm{d} {x} \,\textrm{d} t\\&\quad = - \frac{1}{\varepsilon } \int _0^\tau \int _{\Omega } {\varrho _{\varepsilon }\textbf{u}_\varepsilon \cdot \nabla _{x}G } \ \,\textrm{d} {x} \,\textrm{d} t+ \frac{1}{\varepsilon } \int _0^\tau \int _{\Omega } { \left( 1 - \frac{\overline{\varrho }}{\overline{\varrho } + \varepsilon r} \right) \varrho _{\varepsilon }\textbf{u}_\varepsilon \cdot \nabla _{x}G } \ \,\textrm{d} {x} \,\textrm{d} t\\&\quad = - \frac{1}{\varepsilon } \int _0^\tau \int _{\Omega } {\varrho _{\varepsilon }\textbf{u}_\varepsilon \cdot \nabla _{x}G } \ \,\textrm{d} {x} \,\textrm{d} t+ \int _0^\tau \int _{\Omega } { {\frac{r}{\overline{\varrho } + \varepsilon r}} \varrho _{\varepsilon }\textbf{u}_\varepsilon \cdot \nabla _{x}G } \ \,\textrm{d} {x} \,\textrm{d} t\\&\quad = - \frac{1}{\varepsilon } \int _0^\tau \int _{\Omega } {\varrho _{\varepsilon }\textbf{u}_\varepsilon \cdot \nabla _{x}G } \ \,\textrm{d} {x} \,\textrm{d} t+ \int _0^\tau \int _{\Omega } { r \textbf{u}\cdot \nabla _{x}G } \ \,\textrm{d} {x} \,\textrm{d} t+ {\mathcal {O}}(\varepsilon ). \end{aligned}$$

Consequently, abbreviating $s_\varepsilon =s(\varrho _{\varepsilon }, \vartheta _{\varepsilon })$, relation (6.16) can be rewritten as

$$\begin{aligned}&\int _{\Omega } { E_\varepsilon \left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon \Big | \overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}, \textbf{U}\right) (\tau , \cdot ) } \ \,\textrm{d} {x} \nonumber \\&\qquad + \int _0^\tau \int _{\Omega } { \Big ( {\mathbb {S}} (\overline{\vartheta }, {\mathbb {D}}_x\textbf{u}_\varepsilon ) - {\mathbb {S}} (\overline{\vartheta }, {\mathbb {D}}_x\textbf{U}) \Big ) : \Big ( {\mathbb {D}}_x\textbf{u}_\varepsilon - {\mathbb {D}}_x\textbf{U}\Big ) } \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\qquad +\int _0^\tau \int _{\Omega } { \left( \frac{\overline{\vartheta } + \varepsilon {\mathcal {T}}}{\vartheta _{\varepsilon }^2} \right) \frac{\kappa (\vartheta _{\varepsilon }) \nabla _{x}\vartheta _{\varepsilon }\cdot \nabla _{x}\vartheta _{\varepsilon }}{\varepsilon ^2} } \ \,\textrm{d} {x} \,\textrm{d} t- \int _0^\tau \int _{\Omega } { \frac{\kappa (\vartheta _{\varepsilon }) }{\vartheta _{\varepsilon }} \frac{\nabla _{x}\vartheta _{\varepsilon }}{\varepsilon } \cdot \nabla _{x}{\mathcal {T}}} \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\quad \le - \frac{1}{\varepsilon } \int _0^\tau \int _{\Omega } \left( \varrho _{\varepsilon }(s_\varepsilon - s(\overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}})) \partial _t {\mathcal {T}}\right. \nonumber \\&\qquad \qquad \qquad \qquad \qquad \left. + \varrho _{\varepsilon }(s_\varepsilon - s(\overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}})) \textbf{u}_\varepsilon \cdot \nabla _{x}{\mathcal {T}}\right) \ \,\textrm{d} {x}\,\textrm{d} t\nonumber \\&\qquad - \int _0^\tau \int _{\Omega } { \frac{\varrho _{\varepsilon }}{\varepsilon } \nabla _{x}G \cdot \textbf{U}} \ \,\textrm{d} {x} \,\textrm{d} t+ \int _0^\tau \int _{\Omega } { r \nabla _{x}G \cdot \textbf{U}} \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\qquad + C \int _0^\tau \int _{\Omega } {E_\varepsilon \left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon \Big | \overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}, \textbf{U}\right) } \ \,\textrm{d} {x} \,\textrm{d} t+ {\mathcal {O}}(\varepsilon ). \end{aligned}$$

(6.19)

Step 4: In view of solenoidality $\textrm{div}_{x}\textbf{U}= 0$, we have

$$\begin{aligned}&- \int _0^\tau \int _{\Omega } { \frac{\varrho _{\varepsilon }}{\varepsilon } \nabla _{x}G \cdot \textbf{U}} \ \,\textrm{d} {x} \,\textrm{d} t+ \int _0^\tau \int _{\Omega } { r \nabla _{x}G \cdot \textbf{U}} \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\quad = - \int _0^\tau \int _{\Omega } {\frac{\varrho _{\varepsilon }- (\overline{\varrho } + \varepsilon r)}{\varepsilon } \textbf{U}\cdot \nabla _{x}G } \ \,\textrm{d} {x} \,\textrm{d} t. \end{aligned}$$

(6.20)

Consequently, we may use the bounds (5.20), (5.21) along with the convergences established in (6.1)–(6.6) to rewrite (6.19) in the form

$$\begin{aligned}&\int _{\Omega } { E_\varepsilon \left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon \Big | \overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}, \textbf{U}\right) (\tau , \cdot )} \ \,\textrm{d} {x} \nonumber \\&\qquad + \int _0^\tau \int _{\Omega } { \Big ( {\mathbb {S}} (\overline{\vartheta }, {\mathbb {D}}_x\textbf{u}_\varepsilon ) - {\mathbb {S}} (\overline{\vartheta }, {\mathbb {D}}_x\textbf{U}) \Big ) : \Big ( {\mathbb {D}}_x\textbf{u}_\varepsilon - {\mathbb {D}}_x\textbf{U}\Big ) } \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\qquad +\int _0^\tau \int _{\Omega } { \left( \frac{\overline{\vartheta } + \varepsilon {\mathcal {T}}}{\vartheta _{\varepsilon }^2} \right) \frac{\kappa (\vartheta _{\varepsilon }) \nabla _{x}\vartheta _{\varepsilon }\cdot \nabla _{x}\vartheta _{\varepsilon }}{\varepsilon ^2} } \ \,\textrm{d} {x} \,\textrm{d} t- \int _0^\tau \int _{\Omega } { \frac{\kappa (\overline{\vartheta }) }{\overline{\vartheta }} \nabla _{x}\mathfrak {T} \cdot \nabla _{x}{{\mathcal {T}}} } \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\quad \le - \frac{1}{\varepsilon } \int _0^\tau \int _{\Omega } { \varrho _{\varepsilon }\Big [ (s(\varrho _{\varepsilon }, \vartheta _{\varepsilon }) - s( \overline{\varrho } + \varepsilon r,\overline{\vartheta } + \varepsilon {\mathcal {T}}) \Big ] \partial _t {\mathcal {T}}} \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\qquad - \frac{1}{\varepsilon } \int _0^\tau \int _{\Omega } { \varrho _{\varepsilon }\Big [ s(\varrho _{\varepsilon }, \vartheta _{\varepsilon }) - s( \overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}) \Big ] \textbf{u}_\varepsilon \cdot \nabla _{x}{\mathcal {T}}} \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\qquad + \frac{1}{\varepsilon } \int _0^\tau \int _{\Omega } { \varrho _{\varepsilon }(s(\varrho _{\varepsilon }, \vartheta _{\varepsilon }) - s(\overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}})) ( \textbf{U}- \textbf{u}_\varepsilon ) \cdot \nabla _{x}{\mathcal {T}}} \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\qquad + \int _0^\tau \int _{\Omega } { ( r - \mathfrak {R} ) \nabla _{x}G \cdot \textbf{U}} \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\qquad + C \int _0^\tau \int _{\Omega } {E_\varepsilon \left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon \Big | \overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}, \textbf{U}\right) } \ \,\textrm{d} {x} \,\textrm{d} t+ {\mathcal {O}}(\varepsilon ). \end{aligned}$$

(6.21)

Step 5: Now we use the fact that ${\mathcal {T}}$ solves the modified heat equation (6.10), specifically,

(6.22)

Thus we may perform the limit in several integrals in (6.21) obtaining

$$\begin{aligned}&\int _{\Omega } { E_\varepsilon \left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon \Big | \overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}, \textbf{U}\right) (\tau , \cdot ) } \ \,\textrm{d} {x} \nonumber \\&\qquad + \int _0^\tau \int _{\Omega } { \Big ( {\mathbb {S}} (\overline{\vartheta }, {\mathbb {D}}_x\textbf{u}_\varepsilon ) - {\mathbb {S}} (\overline{\vartheta }, {\mathbb {D}}_x\textbf{U}) \Big ) : \Big ( {\mathbb {D}}_x\textbf{u}_\varepsilon - {\mathbb {D}}_x\textbf{U}\Big ) } \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\qquad +\int _0^\tau \int _{\Omega } { \left( \frac{\overline{\vartheta } + \varepsilon {\mathcal {T}}}{\vartheta _{\varepsilon }^2} \right) \frac{\kappa (\vartheta _{\varepsilon }) \nabla _{x}\vartheta _{\varepsilon }\cdot \nabla _{x}\vartheta _{\varepsilon }}{\varepsilon ^2} } \ \,\textrm{d} {x} \,\textrm{d} t- \int _0^\tau \int _{\Omega } { \frac{\kappa (\overline{\vartheta }) }{\overline{\vartheta }} \nabla _{x}\mathfrak {T} \cdot \nabla _{x}{\mathcal {T}}} \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\quad \le - \int _0^\tau \int _{\Omega } { \overline{\varrho } \left( \frac{\partial s(\overline{\varrho }, \overline{\vartheta })}{\partial \varrho }(\mathfrak {R} - r) + \frac{\partial s(\overline{\varrho }, \overline{\vartheta })}{\partial \vartheta }(\mathfrak {T} - {\mathcal {T}}) \right) {\frac{\overline{\vartheta } \alpha (\overline{\varrho }, \overline{\vartheta } )}{ c_p (\overline{\varrho }, \overline{\vartheta } )} } \nabla _{x}G \cdot \textbf{U}} \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\qquad - { \int _0^\tau \int _{\Omega } { \overline{\varrho } \left( \frac{\partial s(\overline{\varrho }, \overline{\vartheta })}{\partial \varrho }(\mathfrak {R} - r) + \frac{\partial s(\overline{\varrho }, \overline{\vartheta })}{\partial \vartheta }(\mathfrak {T} - {\mathcal {T}}) \right) \frac{\kappa (\overline{\vartheta })}{\overline{\varrho } c_p (\overline{\varrho }, \overline{\vartheta } )} \Delta _x{\mathcal {T}}} \ \,\textrm{d} {x} \,\textrm{d} t} \nonumber \\&\qquad { - \int _0^\tau \int _{\Omega } { \overline{\varrho } \left( \frac{\partial s(\overline{\varrho }, \overline{\vartheta })}{\partial \varrho }(\mathfrak {R} - r) + \frac{\partial s(\overline{\varrho }, \overline{\vartheta })}{\partial \vartheta }(\mathfrak {T} - {\mathcal {T}}) \right) \frac{1}{\overline{\varrho } c_p (\overline{\varrho }, \overline{\vartheta } )} \Lambda (t) } \ \,\textrm{d} {x} \,\textrm{d} t} \nonumber \\&\qquad + \int _0^\tau \int _{\Omega } { ( r - \mathfrak {R} ) \nabla _{x}G \cdot {\textbf{U}} } \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\qquad + C \int _0^\tau \int _{\Omega } {E_\varepsilon \left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon \Big | \overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}, \textbf{U}\right) } \ \,\textrm{d} {x} \,\textrm{d} t+ {\mathcal {O}}(\varepsilon ), \end{aligned}$$

(6.23)

where we have used

$$\begin{aligned}&\frac{1}{\varepsilon } \int _0^\tau \int _{\Omega } { \varrho _{\varepsilon }(s(\varrho _\varepsilon ,\vartheta _\varepsilon ) - s(\overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}})) ( \textbf{U}- \textbf{u}_\varepsilon ) \cdot \nabla _{x}{{\mathcal {T}}} } \ \,\textrm{d} {x} \,\textrm{d} t\\&\quad \lesssim \int _0^\tau \int _{\Omega } {E_\varepsilon \left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon \Big | \overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}, \textbf{U}\right) } \ \,\textrm{d} {x} \,\textrm{d} t. \end{aligned}$$

Now, we use

$$\begin{aligned} \int _{\Omega } { (r - \mathfrak {R}) } \ \,\textrm{d} {x} = 0 \end{aligned}$$

to rewrite the third integral on the right-hand side of (6.23) as

$$\begin{aligned}&- \int _0^\tau \int _{\Omega } { \overline{\varrho } \left( \frac{\partial s(\overline{\varrho }, \overline{\vartheta })}{\partial \varrho }(\mathfrak {R} - r) + \frac{\partial s(\overline{\varrho }, \overline{\vartheta })}{\partial \vartheta }(\mathfrak {T} - {\mathcal {T}}) \right) \frac{1}{\overline{\varrho } c_p (\overline{\varrho }, \overline{\vartheta } )} \Lambda (t)} \ \,\text {d} {x} \,\text {d} t\nonumber \\&\quad =- \int _0^\tau \int _{\Omega } \frac{\partial s(\overline{\varrho }, \overline{\vartheta })}{\partial \vartheta } \left[ \frac{\partial p(\overline{\varrho }, \overline{\vartheta } )}{\partial \varrho } \left( \frac{\partial p(\overline{\varrho }, \overline{\vartheta } )}{\partial \vartheta } \right) ^{-1}(\mathfrak {R} - r) + (\mathfrak {T} - {\mathcal {T}}) \right] \frac{1}{ c_p (\overline{\varrho }, \overline{\vartheta } )} \Lambda (t)\ \,\text {d} {x}\,\text {d} t. \end{aligned}$$

(6.24)

In view of the Boussinesq relations (6.8), (6.11), the expression

$$\begin{aligned} \left[ \frac{\partial p(\overline{\varrho }, \overline{\vartheta } )}{\partial \varrho } \left( \frac{\partial p(\overline{\varrho }, \overline{\vartheta } )}{\partial \vartheta } \right) ^{-1}(\mathfrak {R} - r) + (\mathfrak {T} - {\mathcal {T}}) \right] \end{aligned}$$

is spatially homogeneous, meaning it depends on t only.

Similarly, we can rewrite the second integral on the right-hand side of (6.23) as

$$\begin{aligned}&-\int _0^\tau \int _{\Omega } { \overline{\varrho } \left( \frac{\partial s(\overline{\varrho }, \overline{\vartheta })}{\partial \varrho }(\mathfrak {R} - r) + \frac{\partial s(\overline{\varrho }, \overline{\vartheta })}{\partial \vartheta }(\mathfrak {T} - {\mathcal {T}}) \right) \frac{\kappa (\overline{\vartheta })}{\overline{\varrho } c_p (\overline{\varrho }, \overline{\vartheta } )} \Delta _x{\mathcal {T}}} \ \,\text {d} {x} \,\text {d} t\nonumber \\ {}&\quad = - \int _0^\tau \int _{\Omega } { \frac{\partial s(\overline{\varrho }, \overline{\vartheta })}{\partial \varrho }\left[ (\mathfrak {R} - r) + \frac{\partial p(\overline{\varrho }, \overline{\vartheta }) }{\partial \vartheta } \left( \frac{\partial p(\overline{\varrho }, \overline{\vartheta }) }{\partial \varrho } \right) ^{-1} (\mathfrak {T} - {\mathcal {T}}) \right] \frac{\kappa (\overline{\vartheta })}{ c_p (\overline{\varrho }, \overline{\vartheta } )} \Delta _x{\mathcal {T}}} \ \,\text {d} {x} \,\text {d} t\nonumber \\ {}&\qquad + \int _0^\tau \int _{\Omega } \left( \frac{\partial s(\overline{\varrho }, \overline{\vartheta })}{\partial \varrho }\frac{\partial p(\overline{\varrho }, \overline{\vartheta }) }{\partial \vartheta } \left( \frac{\partial p(\overline{\varrho }, \overline{\vartheta }) }{\partial \varrho } \right) ^{-1} (\mathfrak {T} - {\mathcal {T}}) - \frac{\partial s(\overline{\varrho }, \overline{\vartheta })}{\partial \vartheta } (\mathfrak {T} - {\mathcal {T}})\right) \frac{\kappa (\overline{\vartheta })}{ c_p (\overline{\varrho }, \overline{\vartheta } )} \Delta _x{\mathcal {T}}\ \,\text {d} {x}\,\text {d} t. \end{aligned}$$

(6.25)

Similarly to the above, the quantity

$$\begin{aligned} \left[ (\mathfrak {R} - r) + \frac{\partial p(\overline{\varrho }, \overline{\vartheta }) }{\partial \vartheta } \left( \frac{\partial p(\overline{\varrho }, \overline{\vartheta }) }{\partial \varrho } \right) ^{-1} (\mathfrak {T} - {\mathcal {T}}) \right] \end{aligned}$$

is independent of x.

Now, integrating equation (6.22) in x we obtain the identity

$$\begin{aligned} \left[ \left( \overline{\vartheta } \alpha (\overline{\varrho } ,\overline{\vartheta } ) \frac{\partial p(\overline{\varrho }, \overline{\vartheta })}{\partial \vartheta } \right) ^{-1} -\frac{1}{\overline{\varrho } c_p (\overline{\varrho }, \overline{\vartheta }) } \right] |\Omega | \Lambda (t) = \int _{\partial \Omega } \frac{\kappa (\overline{\vartheta })}{\overline{\varrho } c_p (\overline{\varrho }, \overline{\vartheta })} \nabla _{x}{\mathcal {T}}\cdot \textbf{n} \ \textrm{d}\sigma _x. \end{aligned}$$

(6.26)

Finally, plugging (6.26) into (6.25) we can compute the sum of (6.24) with the first integral in (6.25) obtaining

$$\begin{aligned}&- \int _{\Omega } { \frac{\partial s(\overline{\varrho }, \overline{\vartheta })}{\partial \vartheta } \left[ \frac{\partial p(\overline{\varrho }, \overline{\vartheta } )}{\partial \varrho } \left( \frac{\partial p(\overline{\varrho }, \overline{\vartheta } )}{\partial \vartheta } \right) ^{-1}(\mathfrak {R} - r) + (\mathfrak {T} - {\mathcal {T}}) \right] \frac{1}{ c_p (\overline{\varrho }, \overline{\vartheta } )} \Lambda (t)} \ \,\textrm{d} {x} \nonumber \\&- \int _{\Omega } \frac{\partial s(\overline{\varrho }, \overline{\vartheta })}{\partial \varrho }\left[ (\mathfrak {R} - r) + \frac{\partial p(\overline{\varrho }, \overline{\vartheta }) }{\partial \vartheta } \left( \frac{\partial p(\overline{\varrho }, \overline{\vartheta }) }{\partial \varrho } \right) ^{-1} (\mathfrak {T} - {\mathcal {T}}) \right] \nonumber \\&\quad \times \left[ \overline{\varrho } \left( \overline{\vartheta } \alpha (\overline{\varrho } ,\overline{\vartheta } ) \frac{\partial p(\overline{\varrho }, \overline{\vartheta })}{\partial \vartheta } \right) ^{-1} -\frac{1}{ c_p (\overline{\varrho }, \overline{\vartheta }) } \right] \Lambda (t) \,\textrm{d} {x}. \end{aligned}$$

(6.27)

Now, in accordance with Gibbs’ relation and the definitions of $\alpha$ and $c_p$ in (1.15),

$$\begin{aligned}&- \frac{\partial s(\overline{\varrho }, \overline{\vartheta })}{\partial \vartheta } \frac{\partial p(\overline{\varrho }, \overline{\vartheta } )}{\partial \varrho } \left( \frac{\partial p(\overline{\varrho }, \overline{\vartheta } )}{\partial \vartheta } \right) ^{-1} \nonumber \\&\qquad - \frac{\partial s(\overline{\varrho }, \overline{\vartheta } )}{\partial \varrho } \left( c_p (\overline{\varrho }, \overline{\vartheta })\overline{\varrho } \left( \overline{\vartheta } \alpha (\overline{\varrho } ,\overline{\vartheta } ) \frac{\partial p(\overline{\varrho }, \overline{\vartheta })}{\partial \vartheta } \right) ^{-1} - 1 \right) \nonumber \\&\quad = - \frac{1}{\overline{\vartheta }}\frac{\partial e(\overline{\varrho }, \overline{\vartheta })}{\partial \vartheta } \frac{\partial p(\overline{\varrho }, \overline{\vartheta } )}{\partial \varrho } \left( \frac{\partial p(\overline{\varrho }, \overline{\vartheta } )}{\partial \vartheta } \right) ^{-1} \nonumber \\&\qquad - \frac{\partial s(\overline{\varrho }, \overline{\vartheta } )}{\partial \varrho } \left[ \overline{\varrho } \left( \frac{\partial e(\overline{\varrho }, \overline{\vartheta }) }{\partial \vartheta } + \frac{1}{\overline{\varrho }} \overline{\vartheta } \alpha (\overline{\varrho }, \overline{\vartheta }) \frac{\partial p (\overline{\varrho }, \overline{\vartheta } )}{\partial \vartheta } \right) \left( \overline{\vartheta } \alpha (\overline{\varrho } ,\overline{\vartheta } ) \frac{\partial p(\overline{\varrho }, \overline{\vartheta })}{\partial \vartheta } \right) ^{-1} - 1 \right] \nonumber \\&\quad =- \frac{1}{\overline{\vartheta }}\frac{\partial e(\overline{\varrho }, \overline{\vartheta })}{\partial \vartheta } \frac{\partial p(\overline{\varrho }, \overline{\vartheta } )}{\partial \varrho } \left( \frac{\partial p(\overline{\varrho }, \overline{\vartheta } )}{\partial \vartheta } \right) ^{-1}\nonumber \\&\qquad - \frac{\partial s(\overline{\varrho }, \overline{\vartheta } )}{\partial \varrho } \left[ \overline{\varrho } \frac{\partial e(\overline{\varrho }, \overline{\vartheta }) }{\partial \vartheta } \left( \overline{\vartheta } \alpha (\overline{\varrho } ,\overline{\vartheta } ) \frac{\partial p(\overline{\varrho }, \overline{\vartheta })}{\partial \vartheta } \right) ^{-1}\right] \nonumber \\&\quad =- \frac{1}{\overline{\vartheta }}\frac{\partial e(\overline{\varrho }, \overline{\vartheta })}{\partial \vartheta } \frac{\partial p(\overline{\varrho }, \overline{\vartheta })}{\partial \varrho } \left( \frac{\partial p(\overline{\varrho }, \overline{\vartheta })}{\partial \vartheta } \right) ^{-1} \nonumber \\&\qquad + \frac{1}{\overline{\varrho }} \frac{\partial p(\overline{\varrho }, \overline{\vartheta } )}{\partial \vartheta } \left[ \frac{\partial e(\overline{\varrho }, \overline{\vartheta }) }{\partial \vartheta } \left( \overline{\vartheta } \alpha (\overline{\varrho } ,\overline{\vartheta } ) \frac{\partial p(\overline{\varrho }, \overline{\vartheta })}{\partial \vartheta } \right) ^{-1}\right] = 0. \end{aligned}$$

(6.28)

Thus, the coefficient multiplying $\mathfrak {R} - r$ vanishes. By the same token, we deduce that the coefficient multiplying $\mathfrak {T} - {\mathcal {T}}$ vanishes.

Next, we handle the second integral in (6.25). Using Gibbs’ relation and the constitutive relations obtained in Sect. 2.2, specifically,

$$\begin{aligned} \frac{\partial s (\overline{\varrho }, \overline{\vartheta } ) }{\partial \varrho } = - \frac{1}{\overline{\varrho }^2} \frac{\partial p(\overline{\varrho }, \overline{\vartheta })}{\partial \vartheta }. \end{aligned}$$

Thus, we get

$$\begin{aligned}&\left[ \frac{\partial s(\overline{\varrho }, \overline{\vartheta })}{\partial \varrho }\frac{\partial p(\overline{\varrho }, \overline{\vartheta }) }{\partial \vartheta } \left( \frac{\partial p(\overline{\varrho }, \overline{\vartheta }) }{\partial \varrho } \right) ^{-1} - \frac{\partial s(\overline{\varrho }, \overline{\vartheta })}{\partial \vartheta } \right] \frac{\kappa (\overline{\vartheta })}{ c_p (\overline{\varrho }, \overline{\vartheta } )}\nonumber \\&\quad = -\left[ \frac{1}{\overline{\varrho }^2} \left( \frac{\partial p(\overline{\varrho }, \overline{\vartheta }) }{\partial \vartheta } \right) ^2 \left( \frac{\partial p(\overline{\varrho }, \overline{\vartheta }) }{\partial \varrho } \right) ^{-1} + \frac{1}{\overline{\vartheta }} \frac{\partial e(\overline{\varrho }, \overline{\vartheta })}{\partial \vartheta } \right] \frac{\kappa (\overline{\vartheta })}{ c_p (\overline{\varrho }, \overline{\vartheta } )} ={-\frac{\kappa (\overline{\vartheta })}{\overline{\vartheta }}.} \end{aligned}$$

(6.29)

Finally, we regroup terms containing $\nabla _{x}G$:

$$\begin{aligned}&-\int _{\Omega } { \overline{\varrho } \left( \frac{\partial s(\overline{\varrho }, \overline{\vartheta })}{\partial \varrho }(\mathfrak {R} - r) + \frac{\partial s(\overline{\varrho }, \overline{\vartheta })}{\partial \vartheta }(\mathfrak {T} - {\mathcal {T}}) \right) {\frac{\overline{\vartheta } \alpha (\overline{\varrho }, \overline{\vartheta } )}{ c_p (\overline{\varrho }, \overline{\vartheta } )} } \nabla _{x}G \cdot \textbf{U}} \ \,\textrm{d} {x} \nonumber \\&\qquad + \int _{\Omega } { ( r - \mathfrak {R} ) \nabla _{x}G \cdot {\textbf{U}} } \ \,\textrm{d} {x} \nonumber \\&\quad =\int _{\Omega } { \overline{\varrho } \left( \frac{\partial s(\overline{\varrho }, \overline{\vartheta })}{\partial \varrho }\nabla _{x}(\mathfrak {R} - r) + \frac{\partial s(\overline{\varrho }, \overline{\vartheta })}{\partial \vartheta }\nabla _{x}(\mathfrak {T} - {\mathcal {T}}) \right) {\frac{\overline{\vartheta } \alpha (\overline{\varrho }, \overline{\vartheta } )}{ c_p (\overline{\varrho }, \overline{\vartheta } )} } G \cdot \textbf{U}} \ \,\textrm{d} {x}\nonumber \\&\qquad - \int _{\Omega } { \nabla _{x}( r - \mathfrak {R} ) G \cdot {\textbf{U}} } \ \,\textrm{d} {x}. \end{aligned}$$

(6.30)

Using Boussinesq relation, we deduce

$$\begin{aligned}&\overline{\varrho } \left( \frac{\partial s(\overline{\varrho }, \overline{\vartheta })}{\partial \varrho }\nabla _{x}(\mathfrak {R} - r) + \frac{\partial s(\overline{\varrho }, \overline{\vartheta })}{\partial \vartheta }\nabla _{x}(\mathfrak {T} - {\mathcal {T}}) \right) {\frac{\overline{\vartheta } \alpha (\overline{\varrho }, \overline{\vartheta } )}{ c_p (\overline{\varrho }, \overline{\vartheta } )} } \nonumber \\ {}&\quad =\overline{\varrho } \left( \frac{\partial s(\overline{\varrho }, \overline{\vartheta })}{\partial \varrho }\nabla _{x}(\mathfrak {R} - r) - \frac{\partial s(\overline{\varrho }, \overline{\vartheta })}{\partial \vartheta } {\frac{\partial p(\overline{\varrho }, \overline{\vartheta }) }{\partial \varrho } \left( \frac{\partial p(\overline{\varrho }, \overline{\vartheta }) }{\partial \vartheta } \right) ^{-1} } \nabla _{x}(\mathfrak {R} - r) \right) {\frac{\overline{\vartheta } \alpha (\overline{\varrho }, \overline{\vartheta } )}{ c_p (\overline{\varrho }, \overline{\vartheta } )} } \nonumber \\ {}&\quad =- \overline{\varrho } \left( \frac{1}{\overline{\varrho }^2} \frac{\partial p (\overline{\varrho }, \overline{\vartheta }) }{\partial \vartheta } \nabla _{x}(\mathfrak {R} - r) +\frac{1}{\overline{\vartheta }} \frac{\partial e(\overline{\varrho }, \overline{\vartheta })}{\partial \vartheta } {\frac{\partial p(\overline{\varrho }, \overline{\vartheta }) }{\partial \varrho } \left( \frac{\partial p(\overline{\varrho }, \overline{\vartheta }) }{\partial \vartheta } \right) ^{-1} } \nabla _{x}(\mathfrak {R} - r) \right) {\frac{\overline{\vartheta } \alpha (\overline{\varrho }, \overline{\vartheta } )}{c_p (\overline{\varrho }, \overline{\vartheta } )} }\nonumber \\ {}&\quad = - \nabla _{x}(\mathfrak {R} - r). \end{aligned}$$

(6.31)

Thus, rearranging terms and using $\mathfrak {T}|_{\partial \Omega }={\mathcal {T}}|_{\partial \Omega }$, (6.23) reduces to the desired inequality

$$\begin{aligned}&\int _{\Omega } { E_\varepsilon \left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon \Big | \overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}, \textbf{U}\right) (\tau , \cdot ) } \ \,\textrm{d} {x} \nonumber \\&\qquad + \int _0^\tau \int _{\Omega } { \Big ( {\mathbb {S}} (\overline{\vartheta }, {\mathbb {D}}_x\textbf{u}_\varepsilon ) - {\mathbb {S}} (\overline{\vartheta }, {\mathbb {D}}_x\textbf{U}) \Big ) : \Big ( {\mathbb {D}}_x\textbf{u}_\varepsilon - {\mathbb {D}}_x\textbf{U}\Big ) } \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\qquad + \int _0^\tau \int _{\Omega } { \frac{\kappa (\overline{\vartheta } ) }{\overline{\vartheta }} \left| \nabla _{x}\left( \frac{\vartheta _{\varepsilon }- \overline{\vartheta }}{\varepsilon } \right) - \nabla _{x}{\mathcal {T}}\right| ^2 } \ \,\textrm{d} {x} \,\textrm{d} t\nonumber \\&\quad \lesssim \int _0^\tau \int _{\Omega } {E_\varepsilon \left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_\varepsilon \Big |\overline{\varrho } + \varepsilon r, \overline{\vartheta } + \varepsilon {\mathcal {T}}, \textbf{U}\right) } \ \,\textrm{d} {x} \,\textrm{d} t+ {\mathcal {O}}(\varepsilon ). \end{aligned}$$

(6.32)

Using Grönwall’s lemma and letting $\varepsilon \rightarrow 0$ we obtain the conclusion claimed in Theorem 4.1.