1 Introduction

Consider the motion of a compressible viscous and heat conducting fluid confined between two parallel plates. For simplicity, we suppose the motion is space–periodic with respect to the horizontal variable. Consequently, the spatial domain \(\Omega \) may be identified with

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

The time evolution of the fluid mass density \(\varrho = \varrho (t,x)\), the absolute temperature \(\vartheta = \vartheta (t,x)\), and the velocity \(\textbf{u}= \textbf{u}(t,x)\) is governed by the Navier–Stokes–Fourier (NSF) system:

$$\begin{aligned} \partial _t \varrho + \textrm{div}_x(\varrho \textbf{u})= & {} 0, \end{aligned}$$
(1.1)
$$\begin{aligned} \partial _t (\varrho \textbf{u}) + \textrm{div}_x(\varrho \textbf{u}\otimes \textbf{u}) + \frac{1}{{\varepsilon }^2} \nabla _xp(\varrho , \vartheta )= & {} \textrm{div}_x\mathbb {S}(\vartheta , \nabla _x\textbf{u}) + \frac{1}{{\varepsilon }^{2}} \varrho \nabla _xG, \end{aligned}$$
(1.2)
$$\begin{aligned} \partial _t (\varrho s(\varrho , \vartheta )) + \textrm{div}_x(\varrho s (\varrho , \vartheta ) \textbf{u}) + \textrm{div}_x\left( \frac{ \textbf{q} (\vartheta , \nabla _x\vartheta ) }{\vartheta } \right)= & {} \frac{1}{\vartheta } \left( {\varepsilon }^2 \mathbb {S} : \nabla _x\textbf{u}- \frac{\textbf{q} (\vartheta , \nabla _x\vartheta ) \cdot \nabla _x\vartheta }{\vartheta } \right) , \end{aligned}$$
(1.3)

supplemented with the Dirichlet boundary conditions

$$\begin{aligned} \textbf{u}|_{\partial \Omega }&= 0, \end{aligned}$$
(1.4)
$$\begin{aligned} \vartheta |_{\partial \Omega }&= \vartheta _B. \end{aligned}$$
(1.5)

The viscous stress tensor 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}{3} \textrm{div}_x\textbf{u}\mathbb {I} \right) + \lambda (\vartheta ) \textrm{div}_x\textbf{u}\mathbb {I}, \end{aligned}$$
(1.6)

and the internal energy flux by Fourier’s law

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

The quantity \(s = s(\varrho , \vartheta )\) in (1.3) is the entropy of the system, related to the pressure \(p = p(\varrho , \vartheta )\) and the internal energy \(e = e(\varrho , \vartheta )\) through Gibbs’ equation

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

The potential G represents the effect of gravitation. The Mach number \(\textrm{Ma} = {\varepsilon }\) and the Froude number \(\textrm{Fr} = {\varepsilon }\) are both proportional to a small parameter. If \({\varepsilon }> 0\) is small, the fluid is almost incompressible and strongly stratified, cf. Klein et al. [11]. Our goal is to identify the limit problem for \({\varepsilon }\rightarrow 0\).

1.1 Asymptotic Limit

In accordance with the scaling of (1.2), (1.3), the zero–th order terms in the asymptotic limit are determined by the stationary (static) problem

$$\begin{aligned} \nabla _xp(\varrho , \vartheta ) = \varrho \nabla _xG. \end{aligned}$$
(1.9)

Applying \(\textbf{curl}\) operator to identity (1.9), we successively deduce

$$\begin{aligned} \nabla _x\varrho \times \nabla _xG = 0, \end{aligned}$$

and

$$\begin{aligned} \frac{\partial p (\varrho , \vartheta )}{\partial \varrho } \nabla _x\varrho + \frac{\partial p (\varrho , \vartheta )}{\partial \vartheta } \nabla _x\vartheta = \varrho \nabla _xG \ \Rightarrow \ \nabla _x\vartheta \times \nabla _xG = 0, \end{aligned}$$

where we have anticipated that the pressure also depends non-trivially on the temperature \(\vartheta \) and is such that \(\frac{\partial p(\varrho , \vartheta )}{\partial \vartheta } \ne 0\). Thus for the static problem (1.9) to be solvable, both \(\nabla _x\varrho \) and \(\nabla _x\vartheta \) must be parallel to \(\nabla _xG\). This fact imposes certain restrictions on the distribution of the boundary temperature \(\vartheta _B\). In particular, the motion in an inclined layer studied by Daniels et al. [5] does not admit any static solution. Accordingly, we focus on the particular case

$$\begin{aligned} G = - g x_3,\ \vartheta _B= {\left\{ \begin{array}{ll} \Theta _\text {up} &{} \text { if }\ x_3 = 1, \\ \Theta _\text {bott} &{} \text { if }\ x_3 = 0, \end{array}\right. } \nonumber \\ \text{ where }\ g> 0, \ \Theta _\text {up}> 0, \ \Theta _\text {bott} > 0 \ \text{ are } \text{ constant. } \end{aligned}$$
(1.10)

Fixing the temperature profile \(\vartheta _B = \Theta (x_3)\) to comply with the boundary conditions (1.10), we may recover \(\varrho = r(x_3)\) as a solution of the ODE

$$\begin{aligned} \frac{\partial p (r, \Theta )}{\partial \varrho } \partial _{x_3} r + \frac{\partial p (r, \Theta )}{\partial \vartheta } \partial _{x_3} \Theta = - r g. \end{aligned}$$
(1.11)

Needless to say, such a problem may admit infinitely many solutions.

To simplify, we focus on the case \(\Theta _\textrm{bott} = \Theta _\textrm{up} > 0\). Accordingly, we consider the reference temperature profile \(\Theta = \Theta _\textrm{bott} = \Theta _\textrm{up}\) - a positive constant. Then it follows from (1.11) that the static density profile \(r = r(x_3)\) must be non–constant as long as \(g \ne 0\). Anticipating the asymptotic limit

$$\begin{aligned} \varrho _{\varepsilon }\rightarrow r,\ \vartheta _{\varepsilon }\rightarrow \Theta ,\ \textbf{u}_{\varepsilon }\rightarrow \textbf{U}\ \hbox {(in some sense)} \end{aligned}$$

we deduce from the equation of continuity (1.1)

$$\begin{aligned} \textrm{div}_x(r \textbf{U}) = 0. \end{aligned}$$
(1.12)

Applying (formally) the same argument to the entropy balance (1.3) we get

$$\begin{aligned} \textrm{div}_x(r s(r, \Theta ) \textbf{U}) = 0. \end{aligned}$$
(1.13)

Equations (1.12), (1.13) are compatible only if

$$\begin{aligned} \nabla _xr \cdot \textbf{U}= 0. \end{aligned}$$

As r depends only on the vertical \(x_3\)-variable, this yields

$$\begin{aligned} U_3 \equiv 0 . \end{aligned}$$
(1.14)

In view of the previous arguments, the limit fluid motion exhibits the “stack of pancakes structure” described in Chapter 6 of Majda’s book [14]. Specifically, \(\textbf{U}= [\textbf{U}_h, 0]\), and

$$\begin{aligned}&\frac{\partial p (r, \Theta )}{\partial \varrho } \partial _{x_3} r = - rg, \end{aligned}$$
(1.15)
$$\begin{aligned}&\textrm{div}_h \textbf{U}_h = 0, \end{aligned}$$
(1.16)
$$\begin{aligned}&r \Big ( \partial _t \textbf{U}_h + \textbf{U}_h \cdot \nabla _h \textbf{U}_h \Big ) + \nabla _h \Pi = \mu (\Theta ) \Delta _h \textbf{U}_h + \mu (\Theta ) \partial ^2_{x_3, x_3} \textbf{U}_h. \end{aligned}$$
(1.17)

Here and hereafter, the subscript h refers to the horizontal variable \(x_h = (x_1,x_2)\), \(\nabla _h =[\partial _{x_1}, \partial _{x_2}]\), \(\textrm{div}_h \textbf{v} = \nabla _h \cdot \textbf{v}\), \(\Delta _h = \textrm{div}_h \nabla _h\). The fluid motion is purely horizontal, the coupling between different layers only through the vertical component of the viscous stress.

To the best of our knowledge, there is no rigorous justification of the system (1.15)–(1.17) available in the literature except the inviscid case discussed in [7]. It is worth noting that a similar problem for the barotropic Navier–Stokes system gives rise to a different limit, namely the so–called anelastic approximation, see Masmoudi [15] or Feireisl et al. [8]. Furthermore, as observed in [3], the related case of a low stratification with \({\textrm{Ma}}={\varepsilon }\) and \({\textrm{Fr}} = \sqrt{{\varepsilon }}\) leads to a limiting system of Oberbeck–Boussinesq type with non-local boundary conditions for the temperature.

1.2 The Strategy of the Convergence Proof

We start with the concept of weak solutions for the NSF system with Dirichlet boundary conditions introduced in [4]. In particular, we recall the ballistic energy and the associated relative energy inequality in Sect. 2. Next, we introduce the concept of strong solutions to Majda’s system in Sect. 3. In Sect. 4, we state our main result.

The strategy is of type “weak” \(\rightarrow \) “strong”, meaning the strong solution of the target system is used as a “test function” in the relative energy inequality associated to the primitive system. 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 temperature profiles from their limit values independent of the scaling parameter \({\varepsilon }\). In Sect. 6, we show convergence to the target system (1.15)–(1.17) anticipating the latter admits a regular solution. This formal argument is made rigorous in Sect. 7, where global existence for Majda’s model is established. The last result may be of independent interest.

2 Weak Solutions to the Primitive NSF System

Our analysis is based on the concept of weak solutions to the NSF system introduced in [4], cf. also [10].

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 ; R^3)), \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)

    where \(\vartheta _B\) is an extension of the boundary data to the whole \(\Omega \).

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

    $$\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 \\ {}&= \int _0^T \int _{\Omega } \left[ \mathbb {S}(\vartheta , \nabla _x\textbf{u}) : \nabla _x\varvec{\varphi }- \frac{1}{{\varepsilon }^{2}} \varrho \nabla _xG \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 ; R^3)\).

  • 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}) : \mathbb {D}_x\textbf{u}- \frac{\textbf{q} (\vartheta , \nabla _x\vartheta ) \cdot \nabla _x\vartheta }{\vartheta } \right] \ \,\textrm{d} {x} \,\textrm{d} t + \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\), where \(\mathbb {D}_x\textbf{u}= \frac{1}{2} ( \nabla _x\textbf{u}+ \nabla _x^t \textbf{u})\) is the symmetric gradient.

  • 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 \\ {}&\quad \quad + \int _0^T \psi \int _{\Omega } \frac{\tilde{\vartheta }}{\vartheta } \left[ {\varepsilon }^2 \mathbb {S}(\vartheta , \nabla _x\textbf{u}): \mathbb {D}_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[ \varrho \textbf{u}\cdot \nabla _xG - \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} \,\textrm{d} t \nonumber \\ {}&\quad \quad + \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 })\),

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

2.1 Relative Energy Inequality

In addition to Gibbs’ equation (1.8), 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 \quad \text{ for } \text{ all }\quad \varrho , \vartheta > 0. \end{aligned}$$
(2.6)

Next, following [4], we introduce the scaled relative energy

$$\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 + \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}$$

Now, the hypothesis of thermodynamic stability (2.6) 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

$$\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 - E_{\varepsilon }(\tilde{\varrho }, \widetilde{S}, \widetilde{\textbf{m}}). \end{aligned}$$

Finally, as observed in [4], 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) \ \,\textrm{d} {x} \right] _{t = 0}^{t = \tau } \nonumber \\&\quad \quad + \int _0^\tau \int _{\Omega } \frac{\tilde{\vartheta }}{\vartheta } \left( \mathbb {S} (\vartheta , \nabla _x\textbf{u}) : \mathbb {D}_x\textbf{u}+ \frac{1}{{\varepsilon }^2} \frac{\kappa (\vartheta ) |\nabla _x\vartheta |^2 }{\vartheta } \right) \ \,\textrm{d} {x} \,\textrm{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) \ \,\textrm{d} {x} \,\textrm{d} t \nonumber \\ {}&\quad \quad - \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 , \nabla _x\textbf{u}) \Big ] : \mathbb {D}_x{\tilde{\textbf{u}}} \ \,\textrm{d} {x} \,\textrm{d} t \nonumber \\ {}&\quad \quad + \int _0^\tau \int _{\Omega } \varrho \left[ \frac{1}{{\varepsilon }^{2}} \nabla _xG - \partial _t {\tilde{\textbf{u}}}- ({\tilde{\textbf{u}}}\cdot \nabla _x) {\tilde{\textbf{u}}}\right] \cdot (\textbf{u}- {\tilde{\textbf{u}}}) \ \,\textrm{d} {x} \,\textrm{d} t \nonumber \\ {}&\quad \quad + \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 _xp(\tilde{\varrho }, \tilde{\vartheta }) \right] \ \,\textrm{d} {x} \,\textrm{d} t \end{aligned}$$
(2.7)

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,\quad \tilde{\vartheta }|_{\partial \Omega } = \vartheta _B, \quad {\tilde{\textbf{u}}}|_{\partial \Omega } = 0. \end{aligned}$$
(2.8)

2.2 Constitutive Relations

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

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

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,\quad a > 0. \end{aligned}$$

Accordingly, the internal energy reads

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

Moreover, using several physical principles it was shown in [9, 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,\quad 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, \quad a > 0. \end{aligned}$$
    (2.10)

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

    $$\begin{aligned} P(0) = 0,\ P'(Z)> 0 \quad \text{ for }\quad Z \ge 0,\quad 0 < \frac{ \frac{5}{3} P(Z) - P'(Z) Z }{Z} \le c \quad \text{ for }\quad 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 ),\quad s_\textrm{m} (\varrho , \vartheta ) = \mathcal {S} \left( \frac{\varrho }{\vartheta ^{\frac{3}{2}} } \right) ,\quad 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 [1, 2], 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 ), \ \hbox {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 {\mathop {\sim }\limits ^{<}}\varrho e(\varrho , \vartheta )&{\mathop {\sim }\limits ^{<}}1+ \varrho ^{\frac{5}{3}} + \vartheta ^4, \end{aligned}$$
(2.17)
$$\begin{aligned} s_\textrm{m}(\varrho , \vartheta )&{\mathop {\sim }\limits ^{<}}\left( 1 + |\log (\varrho )| + [\log (\vartheta )]^+ \right) , \end{aligned}$$
(2.18)

see [9, Chapter 3, Section 3.2].

3 Strong Solutions to Majda’s System

Problem (1.16)–(1.17) shares many common features with the \(2d-\)incompressible Navier–Stokes system solved in the celebrated work by Ladyženskaja [12, 13]. Indeed, we show that problem (1.16)–(1.17), endowed with the boundary conditions

$$\begin{aligned} \textbf{U}_h |_{\partial \Omega } = 0,\ \Omega = \mathbb {T}^2 \times (0,1),\ \mathbb {T}^2 = \left( [-1,1] \Big |_{\{ -1,1 \} } \right) ^{2}, \end{aligned}$$
(3.1)

is globally well posed in the framework of Sobolev spaces \(W^{2,p}\) with \(p > 1\) large enough. We report the following result that may be of independent interest.

Theorem 3.1

(Global existence for Majda’s system) Let \(\Theta > 0\) be given. Suppose that

$$\begin{aligned} r \in C^1([0,1]),\ 0 < \underline{r} \le r (x_3) \quad \text{ for } \text{ all }\quad x_3 \in [0,1]. \end{aligned}$$
(3.2)

Let the initial data \(\textbf{U}_{0,h}\) belong to the class

$$\begin{aligned} \textbf{U}_{0,h} \in W^{3,q} \cap W^{1,q}_0(\Omega ; R^2),\ \textrm{div}_h \textbf{U}_{0,h} = 0 \end{aligned}$$
(3.3)

for all \(1 \le q < \infty \).

Then the system (1.16)–(1.17), with the boundary conditions (3.1) and the initial condition (3.3), admits a strong solution \(\textbf{U}_h\) in \((0,T) \times \Omega \), unique in the class

$$\begin{aligned} \partial _t \textbf{U}_h \in L^p(0,T; L^p(\Omega ; R^2)),\quad (\textbf{U}_h , \nabla _h \textbf{U}_h) \in L^p(0,T; W^{2,p}(\Omega ; R^2)\times W^{2,p}(\Omega ; R^{2\times 2} )) \end{aligned}$$
(3.4)

for any \(1 \le p < \infty \).

Remark 3.2

To avoid any misunderstanding we emphasize that by

$$\begin{aligned} \textbf{U}_{0,h} \in W^{3,q} \cap W^{1,q}_0(\Omega ; R^2),\quad \textrm{div}_h \textbf{U}_{0,h} = 0 \end{aligned}$$

for all \(1 \le q < \infty \) we mean

$$\begin{aligned} \textbf{U}_{0,h} \in \bigcap _{q \ge 1} W^{3,q} \cap W^{1,q}_0(\Omega ; R^2),\quad \textrm{div}_h \textbf{U}_{0,h} = 0. \end{aligned}$$

Similarly,

$$\begin{aligned} \partial _t \textbf{U}_h \in L^p(0,T; L^p(\Omega ; R^2)),\quad (\textbf{U}_h, \nabla _h \textbf{U}_h) \in L^p(0,T; W^{2,p}(\Omega ; R^2)\times W^{2,p}(\Omega ; R^{2\times 2} )) \end{aligned}$$

for all finite \(1 \le p < \infty \) means

$$\begin{aligned} \partial _t \textbf{U}_h \in \bigcap _{p \ge 1} L^p(0,T; L^p(\Omega ; R^2)),\quad (\textbf{U}_h, \nabla _h \textbf{U}_h) \in \bigcap _{p \ge 1} L^p(0,T; W^{2,p}(\Omega ; R^2)\times W^{2,p}(\Omega ; R^{2\times 2} )). \end{aligned}$$

The proof of Theorem 3.1 is postponed to Sect. 7.

4 Main Result

Having collected the necessary preliminary material, we are ready to state our main result.

Theorem 4.1

(Singular limit) Let the thermodynamic functions p, e, and s as well as the transport coefficients \(\mu \), \(\lambda \), and \(\kappa \) comply with the structural hypotheses specified in Sect. 2.2. Let

$$\begin{aligned} G = -g x_3,\quad g> 0,\quad \Theta _\textrm{up} = \Theta _\textrm{bott} = \Theta > 0, \end{aligned}$$
(4.1)

and let

$$\begin{aligned} r \in C^1([0,1]) ,\quad 0 < \underline{r} \le r,\quad \frac{\partial p(r, \Theta ) }{\partial \varrho } \partial _{x_3} r = - r g. \end{aligned}$$
(4.2)

Let \((\varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_{\varepsilon })_{{\varepsilon }> 0}\) be a family of weak solutions of the scaled NSF system in the sense of Definition 2.1 emanating from the initial data

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

where

$$\begin{aligned} \int _{\Omega } E_{\varepsilon }\left( \varrho _{0,{\varepsilon }}, \vartheta _{0,{\varepsilon }}, \textbf{u}_{0, {\varepsilon }}\ \Big | \ r, \Theta , [\textbf{U}_{0,h} , 0 ] \right) \ \,\textrm{d} {x} \rightarrow 0 \quad \text{ as }\quad {\varepsilon }\rightarrow 0, \end{aligned}$$
(4.3)

and \(\textbf{U}_{0,h}\) belongs to the class (3.3).

Then

$$\begin{aligned} \textrm{ess} \sup _{\tau \in (0,T)} \int _{\Omega } E_{\varepsilon }\left( \varrho _{{\varepsilon }}, \vartheta _{{\varepsilon }}, \textbf{u}_{{\varepsilon }}\ \Big | \ r, \Theta , [ \textbf{U}_{h}, 0 ] \right) (\tau , \cdot ) \ \,\textrm{d} {x} \rightarrow 0 \quad \text{ as }\quad {\varepsilon }\rightarrow 0, \end{aligned}$$
(4.4)

where \(\textbf{U}_h\) is the unique solution of Majda’s system, the existence of which is guaranteed by Theorem 3.1.

Hypothesis (4.3) corresponds to well–prepared initial data. In view of the coercivity properties of the relative energy stated in (5.1), (5.2) below, relation (4.4) implies, in particular,

$$\begin{aligned} \varrho _{\varepsilon }&\rightarrow r \text{ in }\ L^\infty (0,T; L^{\frac{5}{3}}(\Omega )), \\ \vartheta _{\varepsilon }&\rightarrow \Theta \text{ in }\ L^\infty (0,T; L^2(\Omega )), \\ \varrho _{\varepsilon }\textbf{u}_{\varepsilon }&\rightarrow r [\textbf{U}_h, 0] \ \ \text{ in }\ L^\infty (0,T; L^{1}(\Omega ; R^3)) \end{aligned}$$

as \({\varepsilon }\rightarrow 0\).

The next two sections are devoted to the proof of Theorem 4.1.

5 Uniform Bounds

In order to perform the singular limit in the NSF system we need the associated sequence of weak solutions \((\varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_{\varepsilon })_{{\varepsilon }> 0}\) to be bounded at least in the energy space. First, we introduce the notation of [9] to distinguish between the “essential” and “residual” range of the thermostatic variables \((\varrho , \vartheta )\). Specifically, given a compact set

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

we introduce

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

As shown in [9, 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 E_{{\varepsilon }} \left( \varrho , \vartheta , {\textbf {u}}\Big | \tilde{\varrho }, \tilde{\vartheta }, {\tilde{{\textbf {u}}}}\right) _\text {ess} \ge C \left( \frac{ |\varrho - \tilde{\varrho }|^2 }{{\varepsilon }^2} + \frac{ |\vartheta - \tilde{\vartheta }|^2 }{{\varepsilon }^2} + |{\textbf {u}}- {\tilde{{\textbf {u}}}}|^2 \right) _\text {ess}, \end{aligned}$$
(5.1)
$$\begin{aligned} E_{{\varepsilon }} \left( \varrho , \vartheta , {\textbf {u}}\Big | \tilde{\varrho }, \tilde{\vartheta }, {\tilde{{\textbf {u}}}}\right)&\ge E_{{\varepsilon }} \left( \varrho , \vartheta , {\textbf {u}}\Big | \tilde{\varrho }, \tilde{\vartheta }, {\tilde{{\textbf {u}}}}\right) _\text {res} \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) _\text {res} \end{aligned}$$
(5.2)

whenever \((\tilde{\varrho }, \tilde{\vartheta }) \in \textrm{int}[K]\), where the constant C depends on K and the distance

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

5.1 Energy Estimates for Ill-Prepared Data

We examine a slightly more general situation than in Theorem 4.1. Let \(\Theta >0\) be constant and r the solution of the static problem

$$\begin{aligned} \frac{ \partial p(r, \Theta )}{\partial \varrho } \partial _{x_3} r = - r g. \end{aligned}$$
(5.3)

Next, we consider a family \((\varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_{\varepsilon })_{{\varepsilon }> 0}\) emanating from ill–prepared data \((\varrho _{0,{\varepsilon }}, \vartheta _{0,{\varepsilon }}, \textbf{u}_{0, {\varepsilon }})_{{\varepsilon }> 0}\),

$$\begin{aligned} \int _{\Omega } E_{\varepsilon }\left( \varrho _{0, {\varepsilon }}, \vartheta _{0,{\varepsilon }}, \textbf{u}_{0, {\varepsilon }} \Big | r , \Theta , 0 \right) \ \,\textrm{d} {x} {\mathop {\sim }\limits ^{<}}1 \ \text{ independently } \text{ of }\ {\varepsilon }\rightarrow 0. \end{aligned}$$
(5.4)

The relative energy inequality (2.7) yields

$$\begin{aligned}&\left[ \int _{\Omega } E_{\varepsilon }\left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_{\varepsilon }\Big | r, \Theta , 0 \right) \ \,\textrm{d} {x} \right] _{t = 0}^{t = \tau } \nonumber \\ {}&\qquad + \int _0^\tau \int _{\Omega } \frac{\Theta }{\vartheta _{\varepsilon }} \left( \mathbb {S} (\vartheta _{\varepsilon }, \nabla _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 }^2} \int _0^\tau \int _{\Omega } \frac{\varrho _{\varepsilon }}{r} \left( r\nabla _xG - \nabla _xp(r, \Theta ) \right) \cdot \textbf{u}_{\varepsilon } \ \,\textrm{d} {x} \,\textrm{d} t . \end{aligned}$$
(5.5)

Moreover, in view of (4.1) and (4.2), we deduce the stationary equation

$$\begin{aligned} \nabla _xp(r, \Theta ) = r \nabla _xG; \end{aligned}$$
(5.6)

hence (5.5) reduces to

$$\begin{aligned}&\left[ \int _{\Omega } E_{\varepsilon }\left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_{\varepsilon }\Big | r, \Theta , 0 \right) \ \,\textrm{d} {x} \right] _{t = 0}^{t = \tau } \nonumber \\ {}&\quad + \int _0^\tau \int _{\Omega } \frac{\Theta }{\vartheta _{\varepsilon }} \left( \mathbb {S} (\vartheta _{\varepsilon }, \nabla _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 \le 0. \end{aligned}$$
(5.7)

5.2 Conclusion, Uniform Bounds for Ill-Prepared Data

In view of the estimates obtained in the previous section, we deduce from (5.7) for ill–prepared initial data satisfying (5.4) 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 | r, \Theta , 0 \right) \ \,\textrm{d} {x} {\mathop {\sim }\limits ^{<}}1, \end{aligned}$$
(5.8)
$$\begin{aligned}&\int _0^T \Vert \textbf{u}_{\varepsilon }\Vert ^2_{W^{1,2}_0 (\Omega ; R^3) } \,\textrm{d} t {\mathop {\sim }\limits ^{<}}1, \end{aligned}$$
(5.9)
$$\begin{aligned}&\frac{1}{{\varepsilon }^2} \int _0^T \left( \Vert \nabla _x\log (\vartheta _{\varepsilon }) \Vert ^2_{L^2(\Omega ; R^3)} + \Vert \nabla _x\vartheta _{\varepsilon }^{\frac{\beta }{2}} \Vert ^2_{L^2(\Omega ; R^3)} \right) {\mathop {\sim }\limits ^{<}}1. \end{aligned}$$
(5.10)

Next, it follows from (5.8) 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} {\mathop {\sim }\limits ^{<}}1. \end{aligned}$$
(5.11)

In addition, we get from (5.8):

$$\begin{aligned}&\textrm{ess} \sup _{t \in (0,T)} \int _{\Omega } \varrho _{\varepsilon }|\textbf{u}_{\varepsilon }|^2 \ \,\textrm{d} {x} {\mathop {\sim }\limits ^{<}}1, \nonumber \\ {}&\textrm{ess} \sup _{t \in (0,T)} \left\| \left[ \frac{\varrho _{\varepsilon }- r}{{\varepsilon }} \right] _\textrm{ess} \right\| _{L^2(\Omega )} {\mathop {\sim }\limits ^{<}}1, \nonumber \\ {}&\textrm{ess} \sup _{t \in (0,T)} \left\| \left[ \frac{\vartheta _{\varepsilon }- \Theta }{{\varepsilon }} \right] _\textrm{ess} \right\| _{L^2(\Omega )} {\mathop {\sim }\limits ^{<}}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 )} {\mathop {\sim }\limits ^{<}}1. \end{aligned}$$
(5.12)

Combining (5.10), (5.11), and (5.12), we conclude

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

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 ; R^3)} \,\textrm{d} t {\mathop {\sim }\limits ^{<}}1 \ \text{ for } \text{ some }\ q > 1. \end{aligned}$$
(5.14)

Indeed we have

$$\begin{aligned} \left| \left[ \frac{\kappa (\vartheta _{\varepsilon }) }{\vartheta _{\varepsilon }} \right] _\textrm{res} \frac{\nabla _x\vartheta _{\varepsilon }}{{\varepsilon }} \right| {\mathop {\sim }\limits ^{<}}\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.13). As for the latter, we deduce from (5.10) that

$$\begin{aligned} \left\| \frac{1}{{\varepsilon }} \nabla _x\vartheta _{\varepsilon }^{\frac{\beta }{2}} \right\| _{L^2((0,T) \times \Omega ; R^3)} {\mathop {\sim }\limits ^{<}}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 )} {\mathop {\sim }\limits ^{<}}1 \ \text{ for } \text{ some }\ r > 2. \end{aligned}$$
(5.15)

To see (5.15) first observe that

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

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

$$\begin{aligned} \left\| \vartheta _{\varepsilon }^{\frac{\beta }{2}} \right\| _{L^2(0,T; L^6(\Omega ))} {\mathop {\sim }\limits ^{<}}1. \end{aligned}$$

Consequently, (5.15) follows by interpolation.

Of course, the above uniform bound remain valid also for the well-prepared initial data considered in Theorem 4.1.

6 Convergence to the Target System

We show convergence to the regular solution \(\textbf{U}_h\) in Majda’s system claimed in Theorem 4.1. To get a lean notation, we will identify the two-dimensional velocity \(\textbf{U}_h\) with its three-dimensional counterpart \([\textbf{U}_h, 0]\). The ansatz \((\tilde{\varrho }, \tilde{\vartheta }, {\tilde{\textbf{u}}}) = (r, \Theta , \textbf{U}_h)\) in the relative energy inequality (2.7) yields

$$\begin{aligned}&\left[ \int _{\Omega } E_{\varepsilon }\left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_{\varepsilon }\Big | r , \Theta , \textbf{U}_h \right) \ \,\textrm{d} {x} \right] _{t = 0}^{t = \tau } \nonumber \\ {}&\qquad + \int _0^\tau \int _{\Omega } \frac{\Theta }{\vartheta _{\varepsilon }} \left( \mathbb {S} (\vartheta _{\varepsilon }, \nabla _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 - \int _0^\tau \int _{\Omega } \Big [ \varrho _{\varepsilon }(\textbf{u}_{\varepsilon }- \textbf{U}_h) \otimes (\textbf{u}_{\varepsilon }- \textbf{U}_h) + \frac{1}{{\varepsilon }^2} p(\varrho _{\varepsilon }, \vartheta _{\varepsilon }) \mathbb {I} - \mathbb {S}(\vartheta _{\varepsilon }, \nabla _x\textbf{u}_{\varepsilon }) \Big ] : \mathbb {D}_x\textbf{U}_h \ \,\textrm{d} {x} \,\textrm{d} t \nonumber \\ {}&\qquad + \int _0^\tau \int _{\Omega } \varrho _{\varepsilon }\left[ \partial _t \textbf{U}_h + (\textbf{U}_h \cdot \nabla _x) \textbf{U}_h \right] \cdot (\textbf{U}_h - \textbf{u}_{\varepsilon }) \ \,\textrm{d} {x} \,\textrm{d} t - \frac{1}{{\varepsilon }^2} \int _0^\tau \int _{\Omega } \varrho _{\varepsilon }\nabla _xG \cdot \textbf{U}_h \ \,\textrm{d} {x} \,\textrm{d} t , \end{aligned}$$
(6.1)

where we have used the stationary equation

$$\begin{aligned} \nabla _xp (r, \Theta ) = r \nabla _xG. \end{aligned}$$

Next, seeing that

$$\begin{aligned} \textrm{div}_x\textbf{U}_h = 0,\quad \nabla _xG \cdot \textbf{U}_h = 0, \end{aligned}$$

we deduce

$$\begin{aligned}&\left[ \int _{\Omega } E_{\varepsilon }\left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_{\varepsilon }\Big | r , \Theta , \textbf{U}_h \right) \ \,\textrm{d} {x} \right] _{t = 0}^{t = \tau } \nonumber \\ {}&\qquad + \int _0^\tau \int _{\Omega } \frac{\Theta }{\vartheta _{\varepsilon }} \left( \mathbb {S} (\vartheta _{\varepsilon }, \nabla _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 - \int _0^\tau \int _{\Omega } \Big [ \varrho _{\varepsilon }(\textbf{u}_{\varepsilon }- \textbf{U}_h) \otimes (\textbf{u}_{\varepsilon }- \textbf{U}_h) - \mathbb {S}(\vartheta _{\varepsilon }, \nabla _x\textbf{u}_{\varepsilon }) \Big ] : \mathbb {D}_x\textbf{U}_h \ \,\textrm{d} {x} \,\textrm{d} t \nonumber \\ {}&\qquad + \int _0^\tau \int _{\Omega } \varrho _{\varepsilon }\left[ \partial _t \textbf{U}_h + (\textbf{U}_h \cdot \nabla _x) \textbf{U}_h \right] \cdot (\textbf{U}_h - \textbf{u}_{\varepsilon }) \ \,\textrm{d} {x} \,\textrm{d} t . \end{aligned}$$
(6.2)

Now, in view of the uniform bounds (5.9), (5.12),

$$\begin{aligned} \int _0^\tau&\int _{\Omega } \varrho _{\varepsilon }\left[ \partial _t \textbf{U}_h + (\textbf{U}_h \cdot \nabla _x) \textbf{U}_h \right] \cdot (\textbf{U}_h - \textbf{u}_{\varepsilon }) \ \,\textrm{d} {x} \,\textrm{d} t \\ {}&= \int _0^\tau \int _{\Omega } r \left[ \partial _t \textbf{U}_h + (\textbf{U}_h \cdot \nabla _x) \textbf{U}_h \right] \cdot (\textbf{U}_h - \textbf{u}_{\varepsilon }) \ \,\textrm{d} {x} \,\textrm{d} t + \mathcal {Q}({\varepsilon }), \end{aligned}$$

where \(\mathcal {Q}({\varepsilon })\) denotes a generic function with the property \(\mathcal {Q}({\varepsilon }) \rightarrow 0\) as \({\varepsilon }\rightarrow 0\).

Next, in view of (5.9), (5.12), we may assume

$$\begin{aligned} \varrho _{\varepsilon }\rightarrow r \ \text{ in }\ L^\infty (0,T; L^{\frac{5}{3}}(\Omega )),\ \textbf{u}_{\varepsilon }\rightarrow \textbf{u} \ \text{ weakly } \text{ in }\ L^2(0,T; W^{1,2}_0(\Omega )), \end{aligned}$$

up to a suitable subsequence, where

$$\begin{aligned} \textrm{div}_x(r \textbf{u}) = 0. \end{aligned}$$
(6.3)

Similarly, using the bounds (5.12), (5.13) we may perform the limit in the entropy inequality (2.4) obtaining

$$\begin{aligned} \textrm{div}_x(r s(r, \Theta ) \textbf{u}) \ge 0. \end{aligned}$$

However, thanks to the no–slip boundary conditions,

$$\begin{aligned} \int _{\Omega } \textrm{div}_x(r s(r, \Theta ) \textbf{u}) \ \,\textrm{d} {x} = 0; \end{aligned}$$

therefore

$$\begin{aligned} \textrm{div}_x(r s(r, \Theta ) \textbf{u}) = 0. \end{aligned}$$
(6.4)

Combining (6.3), (6.4) we may infer

$$\begin{aligned} r \frac{ \partial s (r, \Theta ) }{\partial \varrho } \nabla _xr \cdot \textbf{u}= 0. \end{aligned}$$

As entropy is given by the constitutive equation (2.13), (2.14),

$$\begin{aligned} \frac{ \partial s (r, \Theta ) }{\partial \varrho } < 0, \end{aligned}$$

and we conclude

$$\begin{aligned} u_3 = 0,\ \textrm{div}_h \textbf{u}= 0. \end{aligned}$$
(6.5)

Now,

$$\begin{aligned} \int _0^\tau&\int _{\Omega } \varrho _{\varepsilon }\left[ \partial _t \textbf{U}_h + (\textbf{U}_h \cdot \nabla _x) \textbf{U}_h \right] \cdot (\textbf{U}_h - \textbf{u}_{\varepsilon }) \ \,\textrm{d} {x} \,\textrm{d} t \\ {}&= \int _0^\tau \int _{\Omega } r \left[ \partial _t \textbf{U}_h + (\textbf{U}_h \cdot \nabla _x) \textbf{U}_h \right] \cdot (\textbf{U}_h - \textbf{u}) \ \,\textrm{d} {x} \,\textrm{d} t + \mathcal {Q}({\varepsilon }) . \end{aligned}$$

In addition, since \(\textbf{U}_h\), \(\textbf{u}\) satisfy (1.17), (6.5), respectively, we obtain

$$\begin{aligned} \int _0^\tau&\int _{\Omega } r \left[ \partial _t \textbf{U}_h + (\textbf{U}_h \cdot \nabla _x) \textbf{U}_h \right] \cdot (\textbf{U}_h - \textbf{u}) \ \,\textrm{d} {x} \,\textrm{d} t \nonumber \\ {}&= \int _0^\tau \int _{\Omega } \mu (\Theta ) \left[ \Delta _h \textbf{U}_h + \partial ^2_{x_3, x_3} \textbf{U}_h \right] \cdot (\textbf{U}_h - \textbf{u}) \ \,\textrm{d} {x} \,\textrm{d} t \nonumber \\ {}&= - \int _0^\tau \int _{\Omega } \mathbb {S}(\Theta , \nabla _x\textbf{U}_h ) : \mathbb {D}_x(\textbf{U}_h - \textbf{u}) \ \,\textrm{d} {x} \,\textrm{d} t . \end{aligned}$$
(6.6)

Going back to (6.2), we deduce

$$\begin{aligned}&\left[ \int _{\Omega } E_{\varepsilon }\left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, \textbf{u}_{\varepsilon }\Big | r , \Theta , \textbf{U}_h \right) \ \,\textrm{d} {x} \right] _{t = 0}^{t = \tau } \nonumber \\ {}&\qquad + \int _0^\tau \int _{\Omega } \frac{\Theta }{\vartheta _{\varepsilon }} \left( \mathbb {S} (\vartheta _{\varepsilon }, \nabla _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 - \int _0^\tau \int _{\Omega } \Big [ \varrho _{\varepsilon }(\textbf{u}_{\varepsilon }- \textbf{U}_h) \otimes (\textbf{u}_{\varepsilon }- \textbf{U}_h) - \mathbb {S}(\Theta , \nabla _x\textbf{u}) \Big ] : \mathbb {D}_x\textbf{U}_h \ \,\textrm{d} {x} \,\textrm{d} t \nonumber \\ {}&\qquad - \int _0^\tau \int _{\Omega } \mathbb {S}(\Theta , \nabla _x\textbf{U}_h ) : \mathbb {D}_x(\textbf{U}_h - \textbf{u}) \ \,\textrm{d} {x} \,\textrm{d} t + \mathcal {Q}({\varepsilon }) . \end{aligned}$$
(6.7)

Finally, exploiting weak lower semi–continuity of convex functions, we conclude

$$\begin{aligned}&\left[ \int _{\Omega } E_{\varepsilon }\left( \varrho _{\varepsilon }, \vartheta _{\varepsilon }, {\textbf {u}}_{\varepsilon }\Big | r , \Theta , {\textbf {U}}_h \right) \ \,\text {d} {x} \right] _{t = 0}^{t = \tau } \nonumber \\ {}&\qquad + \int _0^\tau \int _{\Omega } \left( \mathbb {S} (\Theta , \nabla _x{\textbf {u}}) - \mathbb {S} (\Theta , \nabla _x{\textbf {U}}_h) \right) : \left( \mathbb {D}_x{\textbf {u}}- \mathbb {D}_x{\textbf {U}}_h \right) \ \,\text {d} {x} \,\text {d} t \nonumber \\ {}&\quad \le - \int _0^\tau \int _{\Omega } \Big [ \varrho _{\varepsilon }({\textbf {u}}_{\varepsilon }- {\textbf {U}}_h) \otimes ({\textbf {u}}_{\varepsilon }- {\textbf {U}}_h) \Big ] : \mathbb {D}_x{\textbf {U}}_h \ \,\text {d} {x} \,\text {d} t + \mathcal {Q}({\varepsilon }) , \end{aligned}$$
(6.8)

which, applying the standard Grönwall argument, yields the desired convergence as well as \(\textbf{u}= \textbf{U}_h\).

We have proved Theorem 4.1.

7 Global Existence for Majda’s Problem

Our ultimate goal is to show global existence of strong solutions to Majda’s model claimed in Theorem 3.1. To this end, it is more convenient to consider the (horizontal) vorticity formulation of (1.16), (1.17). With a slight abuse of notation in the definition of \(\textbf{U}_h\), this formulation reads

$$\begin{aligned} \partial _t \omega + \textbf{U}_h \cdot \nabla _x\omega&= \nu \Delta _x\omega , \end{aligned}$$
(7.1)
$$\begin{aligned} \textbf{U}_h&= \left[ \nabla ^{\perp }_h \Delta _h^{-1}[ \omega ], 0 \right] , \end{aligned}$$
(7.2)
$$\begin{aligned} \nu&= \nu (x_3), \end{aligned}$$
(7.3)

with the boundary conditions

$$\begin{aligned} \omega |_{\partial \Omega } = 0, \end{aligned}$$
(7.4)

and the initial condition

$$\begin{aligned} \omega (0, \cdot ) = \omega _0. \end{aligned}$$
(7.5)

Here, \(\nu = \frac{\mu (\Theta )}{r}\), and

$$\begin{aligned} \omega = \textrm{curl}_h \textbf{U}_h, \ \textrm{curl}_h [\textbf{v}] = \partial _{x_1} v_2 - \partial _{x_2} v_1. \end{aligned}$$
(7.6)

For given \(\omega \), the velocity field \(\textbf{U}_h\) can be recovered via Biot-Savart law:

$$\begin{aligned} \textbf{U}_h = \left[ \nabla ^{\perp }_h \Delta _h^{-1}[ \omega ], 0 \right] ,\ \nabla ^\perp _h = [ -\partial _{x_2},\partial _{x_1} ]. \end{aligned}$$
(7.7)

Remark 7.1

Strictly speaking, the velocity \(\textbf{U}_h\) is determined by (7.7) up to its horizontal average

$$\begin{aligned} \overline{\textbf{U}}_h = \int _{\mathbb {T}^2} \textbf{U}_h \ \textrm{d}x_h \end{aligned}$$

that can be recovered as the unique solution of the parabolic problem

$$\begin{aligned}&r \partial _t \overline{\textbf{U}}_h = \mu (\Theta ) \partial ^2_{x_3,x_3} \overline{\textbf{U}}_h \ \text{ in }\ (0,T) \times (0,1), \\ {}&\overline{\textbf{U}}_h |_{x_3 = 0,1} = 0, \\ {}&\overline{\textbf{U}}_h(0, \cdot ) = \int _{\mathbb {T}^2} \textbf{U}_{0,h} \ \textrm{d}x_h . \end{aligned}$$

7.1 Construction via a Fixed Point Argument

The desired solution \(\omega \) to (7.1)–(7.5) can be constructed via a simple fixed point argument. Consider the set

$$\begin{aligned} X_M = \left\{ \widetilde{\omega } \in C([0,T] \times \overline{\Omega }) \ \Big | \ {\widetilde{\omega }}|_{\partial \Omega } = 0,\ \widetilde{\omega }(0, \cdot ) = \textrm{curl}_h \textbf{U}_{0, h}\,\ \Vert {\widetilde{\omega }} \Vert _{ C([0,T] \times \overline{\Omega }) } \le M \right\} . \end{aligned}$$

As the initial velocity \(\textbf{U}_{0,h}\) belongs to the class (3.3), the set \(X_M\) is a bounded closed convex subset of the Banach space \(C([0,T] \times \overline{\Omega })\). Moreover, \(X_M\) is non-empty as long as M is large enough to accommodate the initial condition.

We define a mapping \(\mathcal {T}[\widetilde{\omega }] = \omega \), where \(\omega \) is the unique solution of the problem

$$\begin{aligned} \partial _t \omega + b_L ( \widetilde{\textbf{U}}_h ) \cdot \nabla _x\omega&= \nu \Delta _x\omega , \end{aligned}$$
(7.8)
$$\begin{aligned} \widetilde{\textbf{U}}_h&= \left[ \nabla ^{\perp }_h \Delta _h^{-1}[ \widetilde{\omega } ], 0 \right] , \end{aligned}$$
(7.9)
$$\begin{aligned} \omega |_{\partial \Omega }&= 0, \end{aligned}$$
(7.10)
$$\begin{aligned} \omega (0, \cdot )&= \textrm{curl}_h \textbf{U}_{0,h} , \end{aligned}$$
(7.11)

for some cut–off function \(b_L\). Specifically,

$$\begin{aligned} b_L(\widetilde{\textbf{U}}_h) = [ b_L( \widetilde{U}_h^1), b_L(\widetilde{U}_h^2), 0], \end{aligned}$$

where

$$\begin{aligned} b_L \in L^\infty (R) \cap C^\infty (R),\ b_L(Z) = Z \ \text{ whenever }\ |Z| \le L. \end{aligned}$$

7.1.1 Maximum Principle

Applying the standard maximum principle, we deduce

$$\begin{aligned} \sup _{t \in [0,T]} \Vert \mathcal {T}[ {\widetilde{\omega }} ](t,\cdot ) \Vert _{C(\overline{\Omega })} = \sup _{t \in [0,T]} \Vert \omega (t,\cdot ) \Vert _{C(\overline{\Omega })} = \Vert \omega (0, \cdot ) \Vert _{C(\overline{\Omega })} {\mathop {\sim }\limits ^{<}}\Vert \textbf{U}_{0,h} \Vert _{W^{2,q} (\Omega ; R^2)} \ \text{ as } \text{ long } \text{ as }\ q > 3. \end{aligned}$$
(7.12)

Note carefully that the bound (7.12) depends solely on the initial data. In particular, it is independent of the specific form of the cut–off function \(b_L\).

7.1.2 Maximal \(L^p-L^q\) Regularity

In view of hypothesis (4.2),

$$\begin{aligned} \nu \in C^1([0,1]),\ 0 < {\underline{\nu }} \le \nu (x_3) \ \text{ for } \text{ any }\ x_3 \in [0,1]. \end{aligned}$$

Consequently, we can apply the maximal \(L^p-L^q\) regularity estimates, see, e.g., Denk, Hieber, and Prüss [6], to obtain

$$\begin{aligned} \Vert \partial _t \omega \Vert _{L^p(0,T; L^q(\Omega ))} + \Vert \omega \Vert _{L^p(0,T; W^{2,q}(\Omega ))}&\le c(p,q) \left( \Vert \omega (0, \cdot ) \Vert _{W^{2,q} \cap W^{1,q}_0(\Omega )} + \Vert b_L (\widetilde{\textbf{U}}_h) \cdot \nabla _x\omega \Vert _{L^p(0,T; L^q(\Omega ))} \right) ,\nonumber \\ 1&< p,q < \infty . \end{aligned}$$
(7.13)

Here

$$\begin{aligned} \Vert \omega (0, \cdot ) \Vert _{W^{2,q} \cap W^{1,q}_0(\Omega )} {\mathop {\sim }\limits ^{<}}\Vert \textbf{U}_{0,h} \Vert _{W^{3,q}(\Omega ; R^2)}, \end{aligned}$$

while, by interpolation and (7.12),

$$\begin{aligned} \Vert b_L (\widetilde{\textbf{U}}_h) \cdot \nabla _x\omega \Vert _{L^q(\Omega )} \le L \Vert \nabla _x\omega \Vert _{L^q(\Omega ; R^3)}&\le L \Vert \omega \Vert _{W^{2,q}(\Omega )}^\lambda \Vert \omega \Vert _{L^q(\Omega )}^{1 - \lambda } \\ {}&\le L c(q) \Vert \textbf{U}_{0,h} \Vert _{W^{3,q}(\Omega ; R^2)}^{1 - \lambda }\Vert \omega \Vert _{W^{2,q}(\Omega )}^\lambda , \ t \in (0,T) \end{aligned}$$

for some \(0< \lambda < 1\). Consequently, it follows from (7.13) and our hypotheses imposed on the initial data that

$$\begin{aligned} \Vert \partial _t \mathcal {T} [ \widetilde{\omega } ] \Vert _{L^p(0,T; L^q(\Omega ))} + \Vert \mathcal {T}[ {\widetilde{\omega }} ] \Vert _{L^p(0,T; W^{2,q}(\Omega ))} \le c\left( p,q, \Vert \textbf{U}_{0,h} \Vert _{W^{3,q}(\Omega ; R^2)} \right) \left( 1 + L \right) \end{aligned}$$
(7.14)

for all finite p, q.

7.2 Fixed Point

It follows from the estimates (7.12), (7.14) that \(\mathcal {T}\) is a compact (continuous) mapping of \(X_M\) into \(X_M\) provided M is large enough, therefore, by means of Tikhonov–Schauder fixed point Theorem, there is a fixed point \(\omega \in X_M\) satisfying

$$\begin{aligned} \partial _t \omega + b_L (\textbf{U}_h ) \cdot \nabla _x\omega&= \nu \Delta _x\omega , \\ \textbf{U}_h&= \left[ \nabla ^{\perp }_h \Delta _h^{-1}[ \omega ], 0 \right] \\ \omega |_{\partial \Omega }&= 0, \\ \omega (0, \cdot )&= \textrm{curl}_h \textbf{U}_{0,h}. \end{aligned}$$

Finally, as \(\textbf{U}_h\) is given by the Biot–Savart law, we get

$$\begin{aligned} \sup _{x_3 \in (0,1)} \Vert \nabla _h \textbf{U}_h \Vert _{L^q(\mathbb {T}^2; R^{2 \times 2})} \le c(q) \Vert \omega (0, \cdot ) \Vert _{L^\infty (\Omega )} \ \text{ uniformly } \text{ for }\ t \in (0,T) \ \text{ for } \text{ any }\ 1< q < \infty , \end{aligned}$$

in particular

$$\begin{aligned} \Vert \textbf{U}_h \Vert _{L^\infty ((0,T) \times \Omega ; R^2)} {\mathop {\sim }\limits ^{<}}\Vert \omega (0, \cdot ) \Vert _{L^\infty (\Omega )} {\mathop {\sim }\limits ^{<}}\Vert \textbf{U}_{0,h} \Vert _{W^{2,q}(\Omega ; R^2)} \ \text{ as } \text{ soon } \text{ as }\ q > 3. \end{aligned}$$

Since this bound is independent of L, we may choose L large enough so that \(b_L (\textbf{U}_h ) = \textbf{U}_h\) to get the desired conclusion

$$\begin{aligned} \partial _t \omega + \textbf{U}_h \cdot \nabla _x\omega&= \nu \Delta _x\omega , \\ \textbf{U}_h&= \left[ \nabla ^{\perp }_h \Delta _h^{-1}[ \omega ], 0 \right] \\ \omega |_{\partial \Omega }&= 0, \\ \omega (0, \cdot )&= \textrm{curl}_h \textbf{U}_{0,h}. \end{aligned}$$

Finally, it is easy to check that the solution is unique in the regularity class (3.4). As a matter of fact, a more general weak–strong uniqueness holds that could be shown adapting the above arguments based on the relative energy inequality.

We have proved Theorem 3.1.