Stationary solutions in thermodynamics of stochastically forced fluids

Abstract

We study the full Navier–Stokes–Fourier system governing the motion of a general viscous, heat-conducting, and compressible fluid subject to stochastic perturbation. The system is supplemented with non-homogeneous Neumann boundary conditions for the temperature and hence energetically open. We show that, in contrast with the energetically closed system, there exists a stationary solution. Our approach is based on new global-in-time estimates which rely on the non-homogeneous boundary conditions combined with estimates for the pressure.

1 Introduction

It is a common belief that the behaviour of turbulent fluid flows can be fully characterized by a steady state of the system (driven by a suitable stochastic forcing to substitute for possible perturbations due to changes in the boundary data), which is approached asymptotically for large times. Mathematically speaking this gives rise to an invariant measure of the underlying system. This is well-understood for the 2D incompressible stochastic Navier–Stokes equations, cf. [10, 12, 16, 19], where uniqueness is well-known.

If uniqueness is not at hand, even the definition of an invariant measure becomes ambiguous, and one rather studies stationary solutions of the dynamics: solutions with a probability law which does not change in time. This law serves as a substitute for an invariant measure. The existence of stationary solutions to the 3D incompressible stochastic Navier–Stokes equations is a nowadays classical result from [11]. More recently a counterpart for the compressible stochastic Navier–Stokes equations has been established in [4]. It is interesting to note that in both cases stationarity provides a certain regularising effect on the solutions (see also [13] in connection with this).

One may think that adding further physical principles such as the possibility of heat transfer completes the picture. The stochastic Navier–Stokes–Fourier equations haven been studied in [1] and the existence of weak martingale solutions has been shown. They describe the motion of a general viscous, heat-conducting, and compressible fluid subject to stochastic perturbation based on the Second Law of Thermodynamics via an entropy balance as in [8] (see also [20] for an alternative approach based on the internal energy balance due to [7]). Supplemented with homogeneous Neumann boundary conditions for the temperature this is an energetically closed system. The mechanical energy which is lost as dissipation is transferred into heat and, different to the incompressible or the isentropic Navier–Stokes equations, weak solutions are known to satisfy an energy equality. The latter one shows that the noise is constantly adding energy to the system such that it can never reach a steady state and, as shown in [1, Section 7], stationary solutions do not exist. Since this is physically not acceptable we are looking for a physical principle which can counteract the energy creation by the noise.

Different to [1] we consider in this paper an energetically open version of the stochastic Navier–Stokes–Fourier equations, where heat can drain through the boundary, see (1.5) below. The time evolution of the fluid in the reference physical domain $Q\subset R^3$ is governed by the following set of equations:

$$\begin{aligned}&\mathrm {d}\varrho +{{\,\mathrm{div}\,}}(\varrho \mathbf{u})\,\mathrm {d}t=0 , \end{aligned}$$

(1.1a)

$$\begin{aligned}&\mathrm {d}(\varrho \mathbf{u})+ \left[ {{\,\mathrm{div}\,}}(\rho \mathbf{u}\otimes \mathbf{u}) + \nabla p(\varrho , \vartheta ) \right] \,\mathrm {d}t= {{\,\mathrm{div}\,}}{\mathbb {S}}(\vartheta , \nabla \mathbf{u}) \,\mathrm {d}t\nonumber \\&\qquad \qquad \qquad \qquad \qquad \qquad \qquad \qquad \qquad \qquad \qquad +\varrho {{\mathbb {F}}}(\varrho ,\vartheta ,\mathbf{u})\,\mathrm {d}W,\,\, \end{aligned}$$

(1.1b)

$$\begin{aligned}&\mathrm {d}(\varrho e(\varrho , \vartheta )) +\bigg [{{\,\mathrm{div}\,}}(\varrho e(\varrho , \vartheta ) \mathbf{u})+{{\,\mathrm{div}\,}}\mathbf{q}(\vartheta , \nabla \vartheta ) \bigg ]\,\mathrm {d}t\nonumber \\&\qquad \qquad \qquad \quad =\bigg [{\mathbb {S}}(\vartheta , \nabla \mathbf{u}) :\nabla \mathbf{u}-p(\varrho , \vartheta ) {{\,\mathrm{div}\,}}\mathbf{u}\bigg ]\,\mathrm {d}t, \end{aligned}$$

(1.1c)

where W is a cylindrical Wiener process and the diffusion coefficient ${\mathbb {F}}$ can be identified with a sequence $({\mathbf {F}}_{k})_{k\ge 1}$ satisfying a suitable Hilbert-Schmidt assumption, see Section 2 for the precise definitions. Here $\varrho$ denotes the density of the fluid, $\vartheta$ the absolute temperature and $\mathbf{u}$ the velocity field. For the viscous stress tensor we suppose Newton’s rheological law

$$\begin{aligned} {\mathbb {S}}={\mathbb {S}}(\vartheta ,\nabla \mathbf{u})=\mu (\vartheta )\Big (\nabla \mathbf{u}+\nabla \mathbf{u}^T-\frac{2}{3}{{\,\mathrm{div}\,}}\mathbf{u}\,{\mathbb {I}}\Big )+\eta (\vartheta ){{\,\mathrm{div}\,}}\mathbf{u}\,{\mathbb {I}}. \end{aligned}$$

(1.2)

The internal energy (heat) flux is determined by Fourier’s law

$$\begin{aligned} \mathbf{q}=\mathbf{q}(\vartheta ,\nabla \vartheta )=-\kappa (\vartheta )\nabla \vartheta =-\nabla {\mathcal {K}}(\vartheta ),\quad {\mathcal {K}}(\vartheta )=\int _0^\vartheta \kappa (z)\,\mathrm {d}z. \end{aligned}$$

(1.3)

The thermodynamic functions p and e are related to the (specific) entropy $s = s(\varrho , \vartheta )$ through Gibbs’ equation

$$\begin{aligned} \vartheta D s (\varrho , \vartheta ) = D e (\varrho , \vartheta ) + p (\varrho , \vartheta ) D \Big ( \frac{1}{\varrho } \Big ) \ \text{ for } \text{ all } \ \varrho , \vartheta > 0, \end{aligned}$$

(1.4)

where D denotes the total derivative with respect to $(\varrho ,\vartheta )$. We supplement (1.1)–(1.4) with the boundary conditions (see also [9])

$$\begin{aligned} \mathbf{u}|_{\partial Q} = 0,\ \mathbf{q}\cdot \mathbf{n}|_{\partial Q} = d(\vartheta )(\vartheta -\Theta _0),\ \text{ and } \text{ fix } \text{ the } \text{ total } \text{ mass }\ \int _Q\varrho \,\mathrm {d}x=M_0, \end{aligned}$$

(1.5)

where $\Theta _0\in L^1(\partial Q)$ is strictly positive, $M_0>0$ and we suppose that there are ${\underline{d}},{\overline{d}}>0$ such that

$$\begin{aligned} {\underline{d}}\vartheta \le d(x,\vartheta )\le {\overline{d}} \vartheta \quad \text {for all}\quad (x,\vartheta )\in \partial Q\times [0,\infty ). \end{aligned}$$

(1.6)

In view of Gibb’s relation (1.4), the internal energy equation (1.1c) can be rewritten in the form of the entropy balance

$$\begin{aligned} \mathrm {d}(\varrho s)+\Big [{{\,\mathrm{div}\,}}(\varrho s\mathbf{u})+{{\,\mathrm{div}\,}}\Big (\frac{\mathbf{q}}{\vartheta }\Big )\Big ]\,\mathrm {d}t=\sigma \,\mathrm {d}t\end{aligned}$$

(1.7)

with the entropy production rate

$$\begin{aligned} \sigma =\frac{1}{\vartheta }\Big ({\mathbb {S}}:\nabla \mathbf{u}-\frac{\mathbf{q}\cdot \nabla \vartheta }{\vartheta }\Big ). \end{aligned}$$

(1.8)

In view of possible singularities, it is convenient to relax the equality sign in (1.8) to the inequality

$$\begin{aligned} \sigma \ge \frac{1}{\vartheta }\Big ({\mathbb {S}}:\nabla \mathbf{u}-\frac{\mathbf{q}\cdot \nabla \vartheta }{\vartheta }\Big ). \end{aligned}$$

(1.9)

The system is augmented by the total energy balance

$$\begin{aligned} \mathrm {d}\int _{Q} \left[ \frac{1}{2} \varrho |\mathbf{u}|^2 + \varrho e \right] \,\mathrm {d}x= & {} \int _{Q} \varrho {\mathbb {F}} \cdot \mathbf{u}\ \mathrm {d}W + \sum _{k\ge 1}\int _{Q} \frac{1}{2} \varrho |\mathbf{F}_k|^2 \,\mathrm {d}x\,\mathrm {d}t\nonumber \\&\quad -\int _{\partial Q}d(\vartheta )(\vartheta -\Theta _0)\,\mathrm {d}{\mathcal {H}}^2\,\mathrm {d}t, \end{aligned}$$

(1.10)

cf. [8, Chapter 2]. In case of a stationary solution applying expectations to (1.10) clearly yields

$$\begin{aligned} \sum _{k\ge 1}{\mathbb {E}}\int _{Q} \frac{1}{2} \varrho |\mathbf{F}_k|^2 \,\mathrm {d}x\,\mathrm {d}t={\mathbb {E}}\int _{\partial Q}d(\vartheta )(\vartheta -\Theta _0)\,\mathrm {d}{\mathcal {H}}^2\,\mathrm {d}t, \end{aligned}$$

meaning energy created by the stochastic forcing can leave through the boundary. The existence theory from [1], which leans on the analysis of the isentropic stochastic Navier–Stokes equations from [5] and the deterministic Navier–Stokes–Fourier equations from [8], can be applied to (1.1)–(1.5) without essential differences. In case of the initial value problem an energy estimate can be derived in terms of the initial data. Looking for stationary solutions, the initial data is not known and one has to use stationarity instead. In [4] stationarity is used in combination with pressure estimates to obtain a corresponding estimate for the isentropic problem. When applying the same strategy to the non-isentropic problem (1.1)–(1.4), supplemented with homogeneous boundary conditions for the temperature flux, the temperature is deemed to grow unboundedly due to the irreversible transfer of the mechanical energy into heat.

Assuming the non-homogeneous boundary conditions (1.5) instead we are able to derive new global-in-time energy estimates, see (4.4). The main task is to control the radiation energy given by $a\vartheta ^4$ without an information on the initial data. In the case of homogeneous boundary conditions one can only obtain informations on the temperature gradient which is not enough to even get estimates for $\vartheta$ in $L^1$. For the non-homogeneous problem we benefit from the boundary term in the energy balance (1.10). A suitable application of Itô’s formula combined with Sobolev’s embedding and an interpolation argument allows to control a higher power of the temperature in terms of the energy, see (4.4). Finally, we derive some pressure estimate by means of the Bogovskii operator in (4.6) and (4.16) to close the argument and to obtain uniform-in-time estimates for the total energy. This leads to our main result which is the existence of stationary martingale solutions to (1.1)–(1.5), see Theorem 2.1 for the precise statement.

In order to make the ideas just explained rigorous one has to regularise the system by adding artificial viscosity to the continuity equation (1.1a) ($\varepsilon$-layer) and add a high power of the pressure in the momentum equation (1.1b) ($\delta$-layer). The resulting system has been solved in [1] by adding three additional layers. The same tedious strategy has been applied in [4] in the construction of stationary solutions to the isentropic system. Here, we follow a different strategy with a much simpler proof. Namely, inspired by the approach due to Itô-Nisio [17] which we recently also applied to the isentropic system with hard sphere pressure [3], we construct stationary solutions directly on the $\varepsilon$-level. The first step is to show uniform-in-time estimates for martingale solutions to the initial value problem. In a second step stationary solutions can be constructed by the Krylov–Bogoliubov method as the narrow limit of time-averages. A striking feature of this approach is that stationary solutions are sitting on the trajectory space and are approached asymptotically in time by any solution starting with bounded initial data of certain moments. With a stationary solution to the approximate system at hand one can prove estimate (4.4) which is uniform in time, $\varepsilon$ and $\delta$. It has to be combined with pressure estimates which differ on both levels, see (4.6) and (4.16), before one can pass to the limit (both limits have to be done independently). The limit passage can be performed as in previous papers and stationarity is preserved in the limit.

2 Mathematical framework and the main result

2.1 Stochastic forcing

The process W is a cylindrical Wiener process on a separable Hilbert space ${\mathfrak {U}}$, that is, $W(t)=\sum _{k\ge 1}\beta _k(t) e_k$ with $(\beta _k)_{k\ge 1}$ being mutually independent real-valued standard Wiener processes relative to $(\mathfrak {F}_t)_{t\ge 0}$. Here $(e_k)_{k\ge 1}$ denotes a complete orthonormal system in $\mathfrak {U}$. In addition, we introduce an auxiliary space $\mathfrak {U}_0\supset \mathfrak {U}$ via

$$\begin{aligned} \mathfrak {U}_0=\bigg \{v=\sum _{k\ge 1}\alpha _k e_k;\;\sum _{k\ge 1}\frac{\alpha _k^2}{k^2}<\infty \bigg \}, \end{aligned}$$

endowed with the norm

$$\begin{aligned} \Vert v\Vert ^2_{\mathfrak {U}_0}=\sum _{k\ge 1}\frac{\alpha _k^2}{k^2},\qquad v=\sum _{k\ge 1}\alpha _k e_k. \end{aligned}$$

Note that the embedding $\mathfrak {U}\hookrightarrow \mathfrak {U}_0$ is Hilbert–Schmidt. Moreover, trajectories of W are ${\mathbb {P}}$-a.s. in $C([0,T];\mathfrak {U}_0)$ (see [6]).

Choosing $\mathfrak {U}= \ell ^2$ we may identify the diffusion coefficients $({\mathbb {F}} e_k )_{k \ge 1}$ with a sequence of real functions $(\mathbf{F}_k)_{k \ge 1}$,

$$\begin{aligned} \varrho {\mathbb {F}}(\varrho , \vartheta , \mathbf{u}) \mathrm {d}W = \sum _{k=1}^\infty \varrho \mathbf{F}_k(x, \varrho , \vartheta , \mathbf{u}) \mathrm {d}\beta _k. \end{aligned}$$

We suppose that $\mathbf{F}_k$ are smooth in their arguments, specifically,

$$\begin{aligned} \mathbf{F}_k \in C^1( \overline{Q} \times [0, \infty )^2 \times R^3; R^3), \end{aligned}$$

where

$$\begin{aligned} \Vert {\mathbf{F}_k} \Vert _{L^\infty } + \Vert {\nabla _{x,\varrho ,\vartheta , \mathbf{u}}} \mathbf{F}_k \Vert _{L^\infty } \le f_k,\ \sum _{k=1}^\infty f_k^2 < \infty . \end{aligned}$$

(2.1)

We easily deduce from (2.1) the following bound

$$\begin{aligned} \Vert \varrho \mathbf{F}_k (\varrho , \vartheta , \mathbf{u}) \Vert _{W^{-k,2}(Q; R^3)} {\mathop {\sim }\limits ^{<}}\Vert \varrho \mathbf{F}_k (\varrho , \vartheta , \mathbf{u}) \Vert _{L^1 (Q; R^3)} {\mathop {\sim }\limits ^{<}}f_k \Vert \varrho \Vert _{L^1(Q)} \end{aligned}$$

whenever $k > \frac{3}{2}$. Accordingly, the stochastic integral

$$\begin{aligned} \int _0^\tau \varrho {\mathbb {F}} \mathrm {d}W = \sum _{k = 1}^\infty \int _0^\tau \varrho \mathbf{F}_k(\varrho , \vartheta , \mathbf{u}) \ \mathrm {d}\beta _k \end{aligned}$$

can be identified with an element of the Banach space space $C([0,T]; W^{-k,2}(Q))$,

$$\begin{aligned} \begin{aligned}&\int _{Q} \left( \int _0^\tau \varrho {\mathbb {F}}(\varrho , \vartheta , \mathbf{u}) \mathrm {d}W \cdot \varvec{\varphi }\right) \ \,\mathrm {d}x \\&\quad = \sum _{k=1}^\infty \int _0^\tau \left( \int _{Q} \varrho \mathbf{F}_k(x,\varrho , \vartheta , \mathbf{u}) \cdot \varphi \ \,\mathrm {d}x \right) \mathrm {d}\beta _k, \ \varvec{\varphi }\in W^{k,2}(Q; R^3),\ k > \frac{3}{2}. \end{aligned} \end{aligned}$$

2.2 Structural and constitutive assumptions

Besides Gibbs’ equation (1.4), we impose several restrictions on the specific shape of the thermodynamic functions $p=p(\varrho ,\vartheta )$, $e=e(\varrho ,\vartheta )$ and $s=s(\varrho ,\vartheta )$. They are borrowed from [8, Chapter 1], to which we refer for the physical background and the relevant discussion.

We consider the pressure p in the form

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

(2.2)

$$\begin{aligned}&e(\varrho , \vartheta ) = e_M(\varrho , \vartheta ) + {a} \frac{\vartheta ^4}{\varrho },\ e_M(\varrho , \vartheta ) = \frac{3}{2} \frac{p_M(\varrho , \vartheta )}{\varrho } = \frac{3}{2} \frac{\vartheta ^{5/2}}{\varrho } P\left( \frac{\varrho }{\vartheta ^{3/2}} \right) , \end{aligned}$$

(2.3)

$$\begin{aligned}&s(\varrho , \vartheta ) = s_M(\varrho , \vartheta ) + \frac{4a }{3} \frac{\vartheta ^3}{\varrho },\ s_M(\varrho , \vartheta ) = S \left( \frac{\varrho }{\vartheta ^{3/2}} \right) , \end{aligned}$$

(2.4)

$$\begin{aligned}&S = S(Z),\ S'(Z)=-\frac{3}{2} \frac{\frac{5}{3} P(Z)-ZP'(Z)}{Z^2}<0,\ \lim _{Z\rightarrow \infty } S(Z)=0, \end{aligned}$$

(2.5)

where

$$\begin{aligned} P \in C^1[0, \infty ) \cap C^2(0, \infty ),\; P(0)= 0,\;P'(Z)>0, \ \text{ for } \text{ all }\ {Z \ge 0}, \end{aligned}$$

(2.6)

$$\begin{aligned} 0<\frac{3}{2} \frac{\frac{5}{3} P(Z)-ZP'(Z)}{Z} <c,\; \text{ for } \text{ all }\ Z > 0, \end{aligned}$$

(2.7)

and

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

(2.8)

As shown in [8, Section 3.2] the assumptions above imply that there is $c>0$ such that

$$\begin{aligned}&c^{-1} \varrho ^{5/3}\le p_M(\varrho ,\vartheta )\le \,c(\varrho ^{5/3}+\varrho \vartheta ), \end{aligned}$$

(2.9)

$$\begin{aligned}&\frac{3p_\infty }{2}\varrho ^{5/3}+ a\vartheta ^4\le \varrho e(\varrho ,\vartheta ), \end{aligned}$$

(2.10)

$$\begin{aligned}&0\le e_M(\varrho ,\vartheta )\le \,{\overline{c}}(\varrho ^{2/3}+\vartheta ), \end{aligned}$$

(2.11)

for all $\vartheta ,\varrho >0$. Moreover, there is $s_\infty >0$ such that

$$\begin{aligned}&0\le s_M(\varrho ,\vartheta )\le \,s_\infty (1+|\log (\varrho )|+ [\log (\vartheta )]^+ ). \end{aligned}$$

(2.12)

Finally, for ${\overline{\vartheta }}>0$ we introduce ballistic free energy given by

$$\begin{aligned} H_{\overline{\vartheta }}(\varrho , \vartheta ) = \varrho \left( e(\varrho , \vartheta ) - \overline{\vartheta } s (\varrho , \vartheta ) \right) , \end{aligned}$$

which satisfies

$$\begin{aligned} -c(\varrho +1)+\frac{1}{4}\big (\varrho e(\varrho ,\vartheta )+{\overline{\vartheta }}|s(\varrho ,\vartheta )|\big )\le H_{\overline{\vartheta }}(\varrho , \vartheta )\le \,c\big (\varrho ^{5/3}+\vartheta ^4+1\big ) \end{aligned}$$

(2.13)

on account of (2.11), (2.12) and [8, Prop. 3.2]. The viscosity coefficients $\mu$, $\eta$ are continuously differentiable functions of the absolute temperature $\vartheta$, more precisely $\mu , \ \lambda \in C^1[0,\infty )$, satisfying

$$\begin{aligned} 0 < {\underline{\mu }}(1 + \vartheta ) \le \mu (\vartheta ) \le \overline{\mu }(1 + \vartheta ), \end{aligned}$$

(2.14)

$$\begin{aligned} \sup _{\vartheta \in [0, \infty )}\big (|\mu '(\vartheta )|+|\lambda '(\vartheta )|\big )\le {\overline{m}}, \end{aligned}$$

(2.15)

$$\begin{aligned} 0 \le \lambda (\vartheta ) \le \overline{\lambda }(1 + \vartheta ). \end{aligned}$$

(2.16)

The heat conductivity coefficient $\kappa \in C^1[0, \infty )$ satisfies

$$\begin{aligned} 0 < {\underline{\kappa }} (1 + \vartheta ^3) \le \kappa ( \vartheta ) \le \overline{\kappa } (1 + \vartheta ^3). \end{aligned}$$

(2.17)

Finally, we introduce certain regularised versions of p, e, s and $\kappa$ for fixed $\delta >0$:

$$\begin{aligned} \begin{aligned} p_\delta (\varrho ,\vartheta )&=p (\varrho ,\vartheta ) + \delta ({\varrho ^2} + \varrho ^\Gamma ), \\ e_{M, \delta }(\varrho , \vartheta )&= e_M(\varrho , \vartheta ) + \delta \vartheta ,\ e_\delta (\varrho , \vartheta ) = e(\varrho , \vartheta ) + \delta \vartheta , \\ s_{M,\delta }(\vartheta ,\varrho )&=s_M(\vartheta ,\varrho )+\delta \log \vartheta ,\quad s_\delta (\varrho , \vartheta ) = s(\varrho , \vartheta ) + \delta \log (\vartheta ),\\ \kappa _\delta (\vartheta )&=\kappa (\vartheta )+\delta \Big (\vartheta ^\Gamma +\frac{1}{\vartheta }\Big ),\ {{\mathcal {K}}}_\delta (\vartheta ) = \int _0^\vartheta \kappa _\delta (z) \ \mathrm{d}z. \end{aligned} \end{aligned}$$

(2.18)

2.3 Martingale and stationary solutions

We start with a rigorous definition of (weak) martingale solution to problem (1.1)–(1.5) as given in [1], where also the existence of a solution to the corresponding initial value problem is proved.

Definition 2.1

(Martingale solution) Let $Q \subset R^3$ be a bounded domain of class $C^{2 + \nu }$, $\nu > 0$. Then

$$\begin{aligned} \big ((\Omega ,\mathfrak {F},(\mathfrak {F}_t),{\mathbb {P}}),\varrho ,\vartheta ,\mathbf{u},W) \end{aligned}$$

is called (weak) martingale solution to problem (1.1)–(1.5) provided the following holds.

• (a) $(\Omega ,\mathfrak {F},(\mathfrak {F}_t),{\mathbb {P}})$ is a stochastic basis with a complete right-continuous filtration;

• (b) W is an $(\mathfrak {F}_t)$-cylindrical Wiener process;

• (c) the random variables

  $$\begin{aligned} \varrho \in L^1_\mathrm{loc}([0,\infty );L^1(Q)),\ \vartheta \in L^1_\mathrm{loc}([0,\infty );L^1(Q)),\ \mathbf{u}\in L^2_\mathrm{loc}([0,\infty ); W^{1,2}_0(Q; R^3)) \end{aligned}$$

  are $(\mathfrak {F}_t)$-progressively measurable¹, $\varrho \ge 0$, $\vartheta > 0$ ${\mathbb {P}}\text{-a.s. }$;

• (d) the equation of continuity

  $$\begin{aligned} \int _0^\infty \int _{Q} \left[ \varrho \partial _t \psi + \varrho \mathbf{u}\cdot \nabla \psi \right] \ \,\mathrm {d}x\,\mathrm {d}t= 0; \end{aligned}$$

  (2.19)

  holds for all $\psi \in C^\infty _c((0,\infty ) \times R^3)$ ${\mathbb {P}}$-a.s.;

• (e) the momentum equation

  $$\begin{aligned}&\int _0^\infty \partial _t \psi \int _{Q} \varrho \mathbf{u}\cdot \varvec{\varphi } \ \,\mathrm {d}x \,\mathrm {d}t\\&\quad +\int _0^\infty \psi \int _{Q} \varrho \mathbf{u}\otimes \mathbf{u}: \nabla \varvec{\varphi } \ \,\mathrm {d}x \,\mathrm {d}t-\int _0^T \psi \int _{Q} {\mathbb {S}}(\vartheta ,\nabla \mathbf{u}) : \nabla \varvec{\varphi } \ \,\mathrm {d}x \,\mathrm {d}t\nonumber \\&\quad +\int _0^\infty \psi \int _{Q} p(\varrho ,\vartheta ) {{\,\mathrm{div}\,}}\varvec{\varphi } \ \,\mathrm {d}x\,\mathrm {d}t+ \int _0^\infty \psi \int _{Q} \varrho {\mathbf{\mathbb {F}}}(\varrho , \vartheta , \mathbf{u}) \cdot \varvec{\varphi } \ \,\mathrm {d}x\,\mathrm {d}W =0;\nonumber \end{aligned}$$

  (2.20)

  holds for all $\psi \in C^\infty _c(0,\infty )$, $\varvec{\varphi }\in C^\infty _c(Q;R^3)$ ${\mathbb {P}}$-a.s.

• (f) the entropy balance

  $$\begin{aligned} \begin{aligned}&- \int _0^\infty \int _{Q} \left[ \varrho s(\varrho ,\vartheta ) \partial _t \psi + \varrho s (\varrho ,\vartheta )\mathbf{u}\cdot \nabla \psi \right] \ \,\mathrm {d}x\mathrm {d}t\\&\ge \int _0^\infty \int _{Q} \frac{1}{\vartheta }\Big [{\mathbb {S}}(\vartheta ,\nabla \mathbf{u}):\nabla \mathbf{u}+\frac{\kappa (\vartheta )}{\vartheta }|\nabla \vartheta |^2\Big ] \psi \,\mathrm {d}x\,\mathrm {d}t \\&- \int _0^\infty \int _{Q} \frac{\kappa (\vartheta )\nabla \vartheta }{\vartheta } \cdot \nabla \psi \ \,\mathrm {d}x\,\mathrm {d}t -\int _0^\infty \int _{\partial Q}\psi \frac{d(\vartheta )}{\vartheta }(\vartheta -\Theta _0)\,\mathrm {d}{\mathcal {H}}^2\,\mathrm {d}t\end{aligned} \end{aligned}$$

  (2.21)

  holds for all $\psi \in C^\infty _c((0,\infty )\times R^3)$, $\psi \ge 0$ ${\mathbb {P}}$-a.s.;

• (g) the total energy balance

  $$\begin{aligned} \begin{aligned} - \int _0^\infty&\partial _t \psi \left( \int _{Q} {\mathcal {E}}(\varrho ,\vartheta , \mathbf{u}) \ \,\mathrm {d}x \right) \ \,\mathrm {d}t= -\int _0^\infty \psi \int _{\partial Q}d(\vartheta -\Theta _0)\,\mathrm {d}{\mathcal {H}}^2\,\mathrm {d}t\\&+\int _0^\infty \psi \int _{Q}\varrho {{\mathbb {F}}}(\varrho ,\vartheta ,\mathbf{u})\cdot \mathbf{u}\,\mathrm {d}x\,\mathrm {d}W\,\mathrm {d}x\\&+ \frac{1}{2} \int _0^\infty \psi \bigg (\int _{Q} \sum _{k \ge 1} \varrho | \mathbf{F}_k (\varrho ,\vartheta , \mathbf{u}) |^2 \ \,\mathrm {d}x \bigg ) \mathrm{d}t \end{aligned} \end{aligned}$$

  (2.22)

  holds for any $\psi \in C^\infty _c(0, \infty )$ ${\mathbb {P}}$-a.s. Here, we abbreviated

  $$\begin{aligned} {\mathcal {E}}(\varrho ,\vartheta ,\mathbf{u})= \frac{1}{2} \varrho | \mathbf{u} |^2 + \varrho e(\varrho ,\vartheta ). \end{aligned}$$

In the following we are going to introduce the concept of stationary martingale solutions. We start with a standard definition of stationarity for stochastic processes with values in Sobolev spaces.

Definition 2.1

(Classical stationarity) Let ${{\tilde{\mathbf{u}}}}=\{{{\tilde{\mathbf{u}}}}(t);t\in [0,\infty )\}$ be an $W^{k,p}(Q)$-valued measurable stochastic process, where $k\in {\mathbb {N}}_0$ and $p\in [1,\infty )$. We say that ${{\tilde{\mathbf{u}}}}$ is stationary on $W^{k,p}(Q)$ provided the joint laws

$$\begin{aligned} {\mathcal {L}}({{\tilde{\mathbf{u}}}}(t_1+\tau ),\dots , {{\tilde{\mathbf{u}}}}(t_n+\tau )),\quad {\mathcal {L}}({{\tilde{\mathbf{u}}}}(t_1),\dots , {{\tilde{\mathbf{u}}}}(t_n)) \end{aligned}$$

on $[W^{k,p}(Q)]^n$ coincide for all $\tau \ge 0$, for all $t_1,\dots ,t_n\in [0,\infty )$.

As can be seen from Definition 2.1, the velocity $\mathbf{u}$ and the temperature $\vartheta$ are not stationary in the sense of Definition 2.1 as they are only equivalence classes in time. Therefore we use the following definition of stationarity which has been introduced in [4], and applies to random variables ranging in the space $L^q_\mathrm{loc}([0, \infty );W^{k,p}(Q))$.

Definition 2.2

(Weak stationarity) Let ${{\tilde{\mathbf{u}}}}$ be an $L^q_\mathrm{loc}([0,\infty );W^{k,p}(Q))$-valued random variable, where $k\in {\mathbb {N}}_0$ and $p,q\in [1,\infty )$. Let ${\mathcal {S}}_\tau$ be the time shift on the space of trajectories given by ${\mathcal {S}}_\tau {{\tilde{\mathbf{u}}}}(t)={{\tilde{\mathbf{u}}}}(t+\tau ).$ We say that ${{\tilde{\mathbf{u}}}}$ is stationary on $L^q_\mathrm{loc}([0,\infty );W^{k,p}(Q))$ provided the laws ${\mathcal {L}}({\mathcal {S}}_\tau {{\tilde{\mathbf{u}}}}),$ ${\mathcal {L}}({{\tilde{\mathbf{u}}}})$ on $L^q_\mathrm{loc}([0,\infty );W^{k,p}(Q))$ coincide for all $\tau \ge 0$.

Definitions 2.1 and 2.2 are equivalent as soon as the stochastic process in question is continuous in time; or alternatively, if it is weakly continuous and satisfies a suitable uniform bound, cf. [4, Lemma A.2 and Corollary A.3]. Furthermore, it can be shown that both notions of stationarity are stable under weak convergence as can be seen from the following two lemmas (the proofs of which can be found in [4, Appendix]).

Lemma 2.1

Let $k\in {\mathbb {N}}_0, p,q\in [1,\infty )$ and let $({{\tilde{\mathbf{u}}}}_m)$ be a sequence of random variables taking values in $L^q_\mathrm{loc}([0,\infty );W^{k,p}(Q)))$. If, for all $m\in {\mathbb {N}}$, ${{\tilde{\mathbf{u}}}}_m$ is stationary on $L^q_\mathrm{loc}([0,\infty );W^{k,p}(Q))$ in the sense of Definition 2.2 and

$$\begin{aligned} {{\tilde{\mathbf{u}}}}_m\rightharpoonup {{\tilde{\mathbf{u}}}}\quad \text {in}\quad L^q_{\mathrm {loc}}([0,\infty );W^{k,p}(Q))\quad {\mathbb {P}}\text {-a.s.,} \end{aligned}$$

then ${{\tilde{\mathbf{u}}}}$ is stationary on $L^q_\mathrm{loc}([0,\infty );W^{k,p}(Q))$.

Lemma 2.2

Let $k\in {\mathbb {N}}_0$, $p\in [1,\infty )$ and let $({{\tilde{\mathbf{u}}}}_m)$ be a sequence of $W^{k,p}(Q)$-valued stochastic processes which are stationary on $W^{k,p}(Q)$ in the sense of Definition 2.1. If for all $T>0$

$$\begin{aligned} \sup _{m\in {\mathbb {N}}} {\mathbb {E}} \left[ \sup _{t\in [0,T]}\Vert {{\tilde{\mathbf{u}}}}_m\Vert _{W^{k,p}(Q)} \right] <\infty \end{aligned}$$

(2.23)

and

$$\begin{aligned} {{\tilde{\mathbf{u}}}}_m\rightarrow {{\tilde{\mathbf{u}}}}\quad \text {in}\quad C_{\mathrm {loc}}([0,\infty );(W^{k,p}(Q),w))\quad {\mathbb {P}}\text {-a.s.,} \end{aligned}$$

then ${{\tilde{\mathbf{u}}}}$ is stationary on $W^{k,p}(Q)$.

In the following we define a stationary martingale solution to (1.1)–(1.5).

Definition 2.3

A weak martingale solution $[\varrho ,\vartheta ,\mathbf{u}, W]$ to (1.1)–(1.5) is called stationary provided the joint law of the time shift $\left[ {\mathcal {S}}_\tau \varrho ,{\mathcal {S}}_\tau \vartheta , {\mathcal {S}}_\tau \mathbf{u}, {\mathcal {S}}_\tau W - W (\tau )\right]$ on

$$\begin{aligned}&L^1_\mathrm{loc}([0, \infty ); L^\gamma ({\mathbb {T}}^3))\times L^1_\mathrm{loc}([0, \infty ); W^{1,2}(Q))\\&\quad \times L^1_\mathrm{loc}([0, \infty ); W^{1,2}(Q; {\mathbb {R}}^3)) \times C([0, \infty ); \mathfrak {U}_0 ) \end{aligned}$$

is independent of $\tau \ge 0$.

We now state our main result concerning the existence of a stationary martingale solution to (1.1)–(1.5).

Theorem 2.1

Let $M_0 \in (0,\infty )$ be given. Suppose that the structural assumptions (2.2)–(2.17) are in force and that the diffusion coecfficient ${\mathbb {F}}$ satisfies (2.1). Then problem (1.1)–(1.5) admits a stationary martingale solution in the sense of Definition 2.3.

The proof of Theorem 2.1 is split into several parts. In the next section we study the approximate system with regularisation parameters $\varepsilon$ and $\delta$. The proof will be completed in Sect. 4 after passing to the limit in $\varepsilon$ and $\delta$.

3 The viscous approximation

In this section we study the viscous approximation to (1.1)–(1.5), where the continuity equation contains an artificial diffusion ($\varepsilon$-layer) and the pressure is stabilised by an artificial high power to the density ($\delta$-layer). In addition to the common terms we add additional stabilising quantities in the continuity equations as in [4], see (3.1) below.

3.1 Martingale solutions

In this subsection we give a precise formulation of the approximated problem. For this purpose we introduce a cut-off function

$$\begin{aligned} \chi \in C^\infty (R), \ \chi (z) = \left\{ \begin{array}{l} 1 \ \text{ for } \ z \le 0, \\ \chi '(z) \le 0 \ \text{ for }\ 0< z < 1, \\ \chi (z) = 0 \ \text{ for }\ z \ge 1. \end{array} \right. \end{aligned}$$

We denote by $M_\varepsilon$ the unique solution to the equation $2\varepsilon z=\chi (z/M_0)$ which obviously satisfies $M_\varepsilon \le M_0$. Finally, the diffusion coefficients are regularised by replacing $\mathbf{\mathbb {F}}$ by $\mathbf{\mathbb {F}}_\varepsilon$,

$$\begin{aligned} {\mathbb {F}}_{\varepsilon } = \left( \mathbf{F}_{k,\varepsilon } \right) _{k \ge 1},\ \mathbf{F}_{k,\varepsilon }(x,\varrho , \vartheta , \mathbf{u}) = \chi \left( \frac{\varepsilon }{\varrho } - 1 \right) \chi \left( |\mathbf{u}|- \frac{1}{\varepsilon } \right) \mathbf{F}_k(x, \varrho , \vartheta , \mathbf{u}). \end{aligned}$$

Let us start with a precise formulation of the problem.

• Regularized equation of continuity.

  $$\begin{aligned} \begin{aligned} \int _{0}^{\infty }&\int _{Q} \left[ \varrho \partial _t \varphi + \varrho \mathbf{u}\cdot \nabla \varphi \right] \ \,\mathrm {d}x \,\mathrm {d}t\\&= \varepsilon \int _0^\infty \int _{Q} \left[ \nabla \varrho \cdot \nabla \varphi - 2 \varrho \varphi \right] \ \,\mathrm {d}x\,\mathrm {d}t- 2\varepsilon \int _0^\infty \int _{Q} M_{\varepsilon } \varphi \ \,\mathrm {d}x \,\mathrm {d}t\end{aligned} \end{aligned}$$

  (3.1)

  for any $\varphi \in C^\infty _c((0, \infty ) \times Q)$ ${\mathbb {P}}$-a.s.

• Regularized momentum equation.

  $$\begin{aligned} \begin{aligned} \int _0^\infty \partial _t \psi&\int _Q \varrho \mathbf{u}\cdot \varphi \,\mathrm {d}x\,\mathrm {d}t+ \int _0^\infty \psi \int _{Q} \varrho \mathbf{u}\otimes \mathbf{u}: \nabla \varphi \ \,\mathrm {d}x \,\mathrm {d}t\\&\quad + \int _0^\infty \psi \int _{Q} p_\delta (\vartheta ,\varrho ) {{\,\mathrm{div}\,}}\varphi \ \,\mathrm {d}x \,\mathrm {d}t\\&\quad - \int _0^\infty \psi \int _{Q} {\mathbb {S}}(\vartheta ,\nabla \mathbf{u}) : \nabla \varphi \ \,\mathrm {d}x \,\mathrm {d}t- \varepsilon \int _0^\infty \psi \int _{Q} \varrho \mathbf{u}\cdot \Delta \varphi \ \,\mathrm {d}x \,\mathrm {d}t\\&\quad - 2 \varepsilon \int _0^\infty \psi \int _{Q} \varrho \mathbf{u}\cdot \varphi \ \,\mathrm {d}x \ \,\mathrm {d}t\\&=- \int _0^\infty \psi \int _{Q} \varrho {\mathbb {F}}_\varepsilon (\varrho ,\vartheta ,\mathbf{u}) \cdot \varphi \ \,\mathrm {d}x \ \mathrm{d}W \end{aligned} \end{aligned}$$

  (3.2)

  for any $\psi \in C^\infty _c((0, \infty ))$, $\varphi \in C^\infty (Q; {\mathbb {R}}^3)$ ${\mathbb {P}}$-a.s.

• Regularized entropy balance.

  $$\begin{aligned} \begin{aligned}&- \int _0^\infty \int _{Q} \left[ \varrho s_\delta (\varrho ,\vartheta ) \partial _t \psi + \varrho s_\delta (\varrho ,\vartheta )\mathbf{u}\cdot \nabla \psi \right] \,\varphi \ \,\mathrm {d}x\mathrm {d}t\\ \ge&\int _0^\infty \int _{Q} \frac{1}{\vartheta }\Big [{\mathbb {S}}(\vartheta ,\nabla \mathbf{u}):\nabla \mathbf{u}+\frac{\kappa _\delta (\vartheta )}{\vartheta }|\nabla \vartheta |^2+\delta \frac{1}{\vartheta ^2}\Big ] \psi \,\mathrm {d}x\,\mathrm {d}t \\&+ \int _0^\infty \psi \int _{Q} \frac{\kappa _\delta (\vartheta )\nabla \vartheta }{\vartheta } \cdot \nabla \varphi \ \,\mathrm {d}x\,\mathrm {d}t \\&-\int _0^\infty \psi \int _{\partial Q}\varphi \frac{d(\vartheta )}{\vartheta }(\vartheta -\Theta _0)\,\mathrm {d}{\mathcal {H}}^2\,\mathrm {d}t\\&-\varepsilon \int _0^\infty \psi \int _Q\left[ \left( \vartheta s_{M,\delta } (\varrho , \vartheta ) - e_{M, \delta } (\varrho , \vartheta ) - \frac{p_M (\varrho , \vartheta ) }{\varrho } \right) \frac{ \nabla \varrho }{\vartheta } \right] \cdot \nabla \varphi \,\mathrm {d}x\,\mathrm {d}t\\&+\int _0^\infty \psi \int _Q\left[ \frac{\varepsilon \delta }{2 \vartheta } ( \beta \varrho ^{\beta - 2} + 2) |\nabla \varrho |^2+\varepsilon \frac{1}{\varrho \vartheta } \frac{\partial p_M}{\partial \varrho } (\varrho , \vartheta ) |\nabla \varrho |^2 - \varepsilon \vartheta ^4\right] \varphi \,\mathrm {d}x\,\mathrm {d}t\\&+\int _0^\infty \psi \int _Q\big (-2\varepsilon \varrho +2\varepsilon M_\varepsilon \big )\frac{1}{\vartheta }\Big ( \vartheta s_{M,\delta }(\varrho , \vartheta ) - e_{M,\delta }(\varrho , \vartheta ) - \frac{p_M(\varrho , \vartheta )}{\varrho } \Big )\,\varphi \,\mathrm {d}x\,\mathrm {d}t\end{aligned} \end{aligned}$$

  (3.3)

  for any $\psi \in C^\infty _c((0, \infty ))$, $\varphi \in C^\infty (Q; {\mathbb {R}}^3)$ ${\mathbb {P}}$-a.s.

• Regularized total energy balance.

  $$\begin{aligned}&\quad - \int _0^\infty \partial _t \psi \bigg (\int _{Q} {{\mathcal {E}}_\delta (\varrho ,\vartheta , \mathbf{u}) } \,\mathrm {d}x\bigg ) \,\mathrm {d}t+ \int _0^\infty {\psi } \int _{Q} \varepsilon \vartheta ^5 \,\mathrm {d}t\nonumber \\&\quad +2 \varepsilon \int _0^T\psi \int _{Q} \left[ \delta \varrho ^2 + \frac{\delta \Gamma }{\Gamma - 1} \varrho ^{\Gamma }+\frac{1}{2}\varrho |\mathbf{u}|^2 \right] \ \,\mathrm {d}x \,\mathrm {d}t\nonumber \\&\quad +\int _0^\infty \psi \int _{\partial Q}d(\vartheta -\Theta _0)\,\mathrm {d}{\mathcal {H}}^2\,\mathrm {d}t\nonumber \\&= \int _0^\infty \int _{Q}\frac{\delta }{\vartheta ^2}\psi \,\mathrm {d}x\,\mathrm {d}t+\int _0^\infty \psi \varepsilon M_\varepsilon \int _{Q} \left( 2\delta \varrho + \frac{\delta \Gamma }{\Gamma - 1} \varrho ^{\Gamma - 1} +\frac{1}{2}|\mathbf{u}|^2\right) \ \,\mathrm {d}x \,\mathrm {d}t\nonumber \\&\quad + \frac{1}{2} \int _0^T \psi \bigg ( \int _{Q} \sum _{k \ge 1} \varrho | \mathbf{F}_{k,\varepsilon } (\varrho ,\vartheta , \mathbf{u}) |^2 \ \,\mathrm {d}x \bigg ) \mathrm{d}t \nonumber \\&\quad +\int _0^\infty \psi \int _{Q}\varrho {\mathbb {F}_\varepsilon }(\varrho ,\vartheta ,\mathbf{u})\cdot \mathbf{u}\,\mathrm {d}W\,\mathrm {d}x\end{aligned}$$

  (3.4)

  holds for any $\psi \in C^\infty _c(0, \infty )$ ${\mathbb {P}}$-a.s., where we have set

  $$\begin{aligned} {\mathcal {E}}_\delta =\frac{1}{2}\varrho |\mathbf{u}|^2+\varrho e_\delta (\varrho ,\vartheta )+\delta \Big (\varrho ^2+\frac{1}{\Gamma -1}\varrho ^\Gamma \Big ). \end{aligned}$$

We have the following result.

Proposition 3.1

Let $\varepsilon ,\delta >0$ be given. Then there exists a weak martingale solution $[\varrho _\varepsilon ,\vartheta _\varepsilon ,\mathbf{u}_\varepsilon ]$ to (3.1)–(3.4). In addition, for $n\in {\mathbb {N}}$ and every $\psi \in C^\infty _c((0,\infty ))$, $\psi \ge 0$, the following generalized energy inequality holds true

$$\begin{aligned} \begin{aligned}&\quad -\int _0^\infty \partial _t\psi \Big [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) \ \,\mathrm {d}x\Big ]^n\,\mathrm {d}t\\&\quad +n\overline{\vartheta } \int _{0}^{\infty }\psi \Big [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Big ]^{n-1}\int _{Q}\sigma _{\varepsilon ,\delta }\,\mathrm {d}x\,\mathrm {d}t\\&\quad + n\int _{0}^{\infty }\psi \Big [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Big ]^{n-1}\int _{Q}\varepsilon \vartheta ^5 \,\mathrm {d}t\\&\quad +n\int _{0}^{\infty }\psi \Big [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Big ]^{n-1}\int _{\partial Q}\frac{d(\vartheta )}{\vartheta }\Big (\vartheta -\overline{\vartheta }\Big )(\vartheta -\Theta _0)\,\mathrm {d}{\mathcal {H}}^2\,\mathrm {d}t\\&\quad + 2\varepsilon n{\overline{\vartheta }}{\mathbb {E}}\int _{0}^{\infty }\psi \Big [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Big ]^{n-1} \int _{Q} \left[ \delta \varrho ^2 + \frac{\delta \Gamma \varrho ^\Gamma }{\Gamma -1} +\frac{1}{2}\varrho |\mathbf{u}|^2\right] \ \,\mathrm {d}x \,\mathrm {d}t\\&\quad +2\varepsilon n \int _{0}^{\infty }\psi \Big [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Big ]^{n-1} \int _{Q}\varrho \Big (\frac{p_M}{\varrho \vartheta }+\frac{e_{M,\delta }}{\vartheta }-s_{M,\delta }\Big )\Big (\varrho ,\vartheta \Big )\,\mathrm {d}x\,\mathrm {d}t\\&{\le }n\varepsilon \int _{0}^{\infty }\psi \Big [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Big ]^{n-1} \int _{Q}\frac{\overline{\vartheta }}{\vartheta ^2}\bigg (e_{M,\delta }(\varrho ,\vartheta )+\varrho \frac{\partial e_M}{\partial \varrho }(\varrho ,\vartheta )\bigg )\nabla \varrho \cdot \nabla \vartheta \,\mathrm {d}x\,\mathrm {d}t\\&\quad +n\int _{0}^{\infty }\psi \Big [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Big ]^{n-1} \int _{Q}\varrho {\mathbf{F}}_\varepsilon (\varrho ,\vartheta ,\mathbf{u})\cdot \mathbf{u}\,\mathrm {d}W\,\mathrm {d}x\\&\quad + \frac{n}{2}\int _{0}^{\infty }\psi \Big [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Big ]^{n-1} \bigg ( \int _{Q} \sum _{k \ge 1} \varrho | \mathbf{F}_{k,\varepsilon } (\varrho ,\vartheta ,\mathbf{u}) |^2 \ \,\mathrm {d}x \bigg ) \mathrm{d}t\\&\quad + n\varepsilon M_\varepsilon \int _{0}^{\infty }\psi \Big [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Big ]^{n-1} \int _{Q} \left( 2\delta \varrho + \frac{\delta \Gamma \varrho ^{\Gamma -1}}{\Gamma -1}+\frac{1}{2}|\mathbf{u}|^2\right) \ \,\mathrm {d}x \,\mathrm {d}t\\&\quad + n{\overline{\vartheta }}\varepsilon M_\varepsilon \int _{0}^{\infty }\psi \Big [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Big ]^{n-1} \int _{Q} \left( \frac{p_M}{\varrho \vartheta }+\frac{e_{M,\delta }}{\vartheta }-s_{M,\delta }\Big )\Big (\varrho ,\vartheta \right) \ \,\mathrm {d}x \,\mathrm {d}t\\&\quad +\frac{n(n-1)}{2}\int _{0}^{\infty }\psi \Big [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Big ]^{n-2} \sum _{k=1}^\infty \int _{\tau _1}^{\tau _2}\bigg ( \int _{Q}\varrho \mathbf{F}_{k,\varepsilon } (\varrho ,\vartheta ,\mathbf{u})\cdot \mathbf{u}\,\mathrm {d}x\bigg )^2\,\mathrm {d}t.\\ \end{aligned} \end{aligned}$$

(3.5)

Here we abbreviated

$$\begin{aligned} \sigma _{\varepsilon ,\delta }&=\frac{1}{\vartheta }\Big [{\mathbb {S}}(\vartheta ,\nabla \mathbf{u}):\nabla \mathbf{u}+\frac{\kappa (\vartheta )}{\vartheta }|\nabla \vartheta |^2+\frac{\delta }{2}\Big (\vartheta ^{\Gamma -1}+\frac{1}{\vartheta ^2}\Big )|\nabla \vartheta |^2+\delta \frac{1}{\vartheta ^2}\Big ]\\&\quad +\frac{\varepsilon \delta }{2\vartheta } {\left( \beta \varrho ^{\Gamma -2} + 2 \right) } |\nabla \varrho |^2+\varepsilon \frac{\partial p_M}{\partial \varrho }(\varrho ,\vartheta )\frac{|\nabla \varrho |^2}{\varrho \vartheta } + \varepsilon \frac{\varrho }{\vartheta } |\nabla \mathbf{u}|^2 , \end{aligned}$$

and $E_{H}^{\delta , \overline{\vartheta }}=\frac{1}{2}\varrho |\mathbf{u}|^2+H_{\overline{\vartheta }}(\varrho , \vartheta )+\delta \big (\varrho ^2+\frac{1}{\Gamma -1}\varrho ^\Gamma \big )$, where

$$\begin{aligned} H_{\delta , \overline{\vartheta }}(\varrho , \vartheta ) = \varrho \left( e_\delta (\varrho , \vartheta ) - \overline{\vartheta } s_\delta (\varrho , \vartheta ) \right) =H_{{\overline{\vartheta }}}(\varrho , \vartheta )+\delta \varrho \vartheta -{\overline{\vartheta }} \varrho \log (\vartheta ) \end{aligned}$$

(3.6)

with $H_{\overline{\vartheta }}(\varrho , \vartheta )$ introduced in (2.2).

Proof

Although there are some differences to system (4.24)–(4.27) from [1] the method still applies (in particular, it is possible to allow an unbounded time interval by working with spaces of the from $L^q_\mathrm{loc}([0,\infty );X)$ and $C_\mathrm{loc}([,0,\infty );X)$ for Banach spaces X) and we obtain the existence of a weak martingale solution to (3.1)–(3.4). We remark, in particular, that the solution in [1] is constructed with respect to some initial law which does not play any role in our analysis. For simplicity we choose

$$\begin{aligned} \varrho _0=1,\quad \vartheta _0=1,\quad \mathbf{u}_0=0, \end{aligned}$$

which satisfies all the required assumptions.

As far as the energy inequality is concerned, the required version can be derived on the basic approximate level (even with equality) and it is preserved in the limit. In fact, one can argue as in [1, Section 4.1] to derive the version for $n=1$, while the case $n\ge 2$ follows easily from Itô’s formula. It is worth to point out that this procedure to test the continuity equation (3.1) with $\frac{1}{2}|\mathbf{u}|^2$ and $2\delta \varrho +\frac{\delta \Gamma }{\Gamma -1}\varrho ^\Gamma$ gives rise to the terms

$$\begin{aligned} 2\varepsilon n{\overline{\vartheta }}{\mathbb {E}}\int _{0}^{\infty }\psi \Big [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Big ]^{n-1}\int _{Q} \left[ \delta \varrho ^2 + \frac{\delta \Gamma }{\Gamma - 1} \varrho ^{\Gamma } +\frac{1}{2}\varrho |\mathbf{u}|^2\right] \ \,\mathrm {d}x \,\mathrm {d}t\end{aligned}$$

and

$$\begin{aligned} n\varepsilon M_\varepsilon \int _{0}^{\infty }\psi \Big [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Big ]^{n-1} \int _{Q} \left( 2\delta \varrho + \frac{\delta \Gamma }{\Gamma - 1} \varrho ^{\Gamma - 1}+\frac{1}{2}|\mathbf{u}|^2\right) \ \,\mathrm {d}x \,\mathrm {d}t\end{aligned}$$

in (3.5) which are new in comparison to [1]. Also, the term

$$\begin{aligned}&n{\overline{\vartheta }}\varepsilon M_\varepsilon \int _{0}^{\infty }\psi \Big [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Big ]^{n-1} \\&\quad \int _{Q} \left( \frac{p_M(\varrho , \vartheta )}{\varrho \vartheta }+\frac{e_{M,\delta }(\varrho , \vartheta )}{\vartheta } - s_{M,\delta }(\varrho , \vartheta ) \right) \ \,\mathrm {d}x \,\mathrm {d}t, \end{aligned}$$

which arises due to the last line in (3.3), does not appear in [1]. Finally, as in (1.10) we have the boundary term due to non-homogeneous boundary conditions being incorporated already in (3.3). $\square$

3.2 Uniform-in-time estimates

The first step is now to derive estimates which are uniform in time.

Proposition 3.2

Let $(\varrho ,\vartheta ,\mathbf{u})$ be a weak martingale solution to (3.1)–(3.4). Assume that

$$\begin{aligned} \mathrm{ess} \limsup _{t \rightarrow 0+} {\mathbb {E}}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (t, \cdot ) \ \,\mathrm {d}x\Bigg ]^n<\infty \end{aligned}$$

(3.7)

for some $n\in {\mathbb {N}}$. Then for any ${\overline{\vartheta }}>0$ and $\varepsilon \le \varepsilon _0$ there is $E_\infty =E_\infty (n,\varepsilon ,\delta ,\vartheta )$ such that

$$\begin{aligned} {\mathbb {E}}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^n&\le E_\infty , \end{aligned}$$

(3.8)

as well as

$$\begin{aligned}&{\mathbb {E}}\int _{0}^{\tau }\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n-1}\bigg (\int _{Q}\sigma _{\varepsilon ,\delta }\,\mathrm {d}x+\int _Q\varepsilon \vartheta ^5\,\mathrm {d}x\bigg ) \,\mathrm {d}t\nonumber \\&\le \,E_\infty (1+\tau ) \end{aligned}$$

(3.9)

for a.a. $\tau > 0$.

Proof

The energy inequality in (3.5) yields for a.a. $0< \tau _1<\tau _2$ (approximating $\psi ={\mathbb {I}}_{(\tau _1,\tau _2)}$ by a sequence of functions $(\psi _m)\subset C^\infty _c(0,\infty )$)

$$\begin{aligned}&\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau _2, \cdot ) \ \,\mathrm {d}x\Bigg ]^n +n\overline{\vartheta } \int _{\tau _1}^{\tau _2}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n-1}\int _{Q}\sigma _{\varepsilon ,\delta }\,\mathrm {d}x\nonumber \\&+ n\int _{\tau _1}^{\tau _2}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n-1}\bigg (\int _{Q}\varepsilon \vartheta ^5\,\mathrm {d}x+\int _{\partial Q}\frac{d(\vartheta )}{\vartheta }\Bigg (\vartheta -\overline{\vartheta }\Bigg )(\vartheta -\Theta _0)\,\mathrm {d}{\mathcal {H}}^2\bigg ) \,\mathrm {d}t\nonumber \\&+ 2\varepsilon n{\overline{\vartheta }}{\mathbb {E}}\int _{\tau _1}^{\tau _2}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n-1} \int _{Q} \left[ \delta \varrho ^2 + \frac{\delta \Gamma }{\Gamma - 1} \varrho ^{\Gamma } +\frac{1}{2}\varrho |\mathbf{u}|^2\right] \ \,\mathrm {d}x \,\mathrm {d}t\nonumber \\&+2\varepsilon n \int _{\tau _1}^{\tau _2} \Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n-1} \int _{Q}\varrho \Bigg (e_{M,\delta }+\varrho \frac{\partial e_{M,\delta }}{\partial \varrho }\Big )\Big (\varrho ,\vartheta \Big )\,\mathrm {d}x\,\mathrm {d}t\nonumber \\&{\le } \,\Bigg [\int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta )(\tau _1,\cdot ) \ \,\mathrm {d}x\Bigg ]^n\nonumber \\&+n\varepsilon \int _{\tau _1}^{\tau _2}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n-1}\int _{Q}\frac{\overline{\vartheta }}{\vartheta ^2}\bigg (e_{M,\delta }+\varrho \frac{\partial e_{M,\delta }}{\partial \varrho })(\varrho ,\vartheta )\bigg ) \nabla \varrho \cdot \nabla \vartheta \,\mathrm {d}x\,\mathrm {d}t\nonumber \\&+n\int _{\tau _1}^{\tau _2}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n-1}\int _{Q}\bigg (\frac{\delta }{\vartheta ^2}+\varepsilon \overline{\vartheta } \vartheta ^4\bigg )\,\mathrm {d}x\,\mathrm {d}t\nonumber \\&+n\int _{\tau _1}^{\tau _2}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n-1} \int _{Q}\varrho {{\mathbb {F}}}_{\varepsilon } (\varrho ,\vartheta ,\mathbf{u})\cdot \mathbf{u}\,\mathrm {d}W\,\mathrm {d}x\nonumber \\&+ \frac{n}{2}\int _{\tau _1}^{\tau _2} \Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n-1}\bigg ( \int _{Q} \sum _{k \ge 1} \varrho | \mathbf{F}_{k,\varepsilon } (\varrho ,\vartheta ,\mathbf{u}) |^2 \ \,\mathrm {d}x \bigg ) \mathrm{d}t \nonumber \\&+ \varepsilon M_\varepsilon n\int _{\tau _1}^{\tau _2} \Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n-1} \int _{Q} \left( 2\delta \varrho + \frac{\delta \Gamma \varrho ^{\Gamma -1}}{\Gamma -1}+\frac{1}{2}|\mathbf{u}|^2\right) \ \,\mathrm {d}x \,\mathrm {d}t\nonumber \\&+ \varepsilon M_\varepsilon n{\overline{\vartheta }}\int _{\tau _1}^{\tau _2} \Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n-1} \int _{Q} \left( \frac{p_M}{\varrho \vartheta }+\frac{e_{M,\delta }}{\vartheta }-s_{M,\delta }\right) \Big (\varrho ,\vartheta \Big ) \ \,\mathrm {d}x \,\mathrm {d}t\nonumber \\&+\frac{n(n-1)}{2}\int _{\tau _1}^{\tau _2} \Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n-2}\sum _{k=1}^\infty \bigg ( \int _{Q}\varrho {\mathbf{F}_{k,\varepsilon }} (\varrho ,\vartheta ,\mathbf{u})\cdot \mathbf{u}\,\mathrm {d}x\bigg )^2\,\mathrm {d}t\nonumber \\&=:(I)+(II)+\dots + (VIII). \end{aligned}$$

(3.10)

Let us first consider the terms on the left-hand side. We have

$$\begin{aligned} n\int _{\tau _1}^{\tau _2}&\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n-1}\int _{\partial Q}\frac{d(\vartheta )}{\vartheta }\Bigg (\vartheta -\overline{\vartheta }\Bigg )(\vartheta -\Theta _0)\,\mathrm {d}{\mathcal {H}}^2\,\mathrm {d}t\\&\ge n\int _{\tau _1}^{\tau _2}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n-1}\int _{\partial Q}\vartheta ^2\,\mathrm {d}{\mathcal {H}}^2\,\mathrm {d}t\\&\quad -cn\int _{\tau _1}^{\tau _2}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n-1}\\&\ge n\int _{\tau _1}^{\tau _2}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n-1}\int _{\partial Q}\vartheta ^2\,\mathrm {d}{\mathcal {H}}^2\,\mathrm {d}t\\ {}&\quad -\kappa n\int _{\tau _1}^{\tau _2}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n}\,\mathrm {d}t-c_\kappa (\tau _2-\tau _1) \end{aligned}$$

for all $\kappa >0$ as well as

$$\begin{aligned} \int _{\tau _1}^{\tau _2}&\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n-1}\int _{Q} \left[ \delta \varrho ^2 + \frac{\delta \Gamma }{\Gamma - 1} \varrho ^{\Gamma } +\frac{1}{2}\varrho |\mathbf{u}|^2\right] \ \,\mathrm {d}x \,\mathrm {d}t\\&\ge {\underline{c}}\int _{\tau _1}^{\tau _2}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n}\,\mathrm {d}t-c(\tau _2-\tau _1)\\ {}&\quad -c\int _{\tau _1}^{\tau _2}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n-1}\int _Q \bigg (\vartheta ^4+{\overline{\vartheta }}\varrho \log (\vartheta )\bigg ) \,\mathrm {d}t\\&\ge {\underline{c}}\int _{\tau _1}^{\tau _2}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n}\,\mathrm {d}t-c_\kappa (\tau _2-\tau _1)\\ {}&\quad -\kappa \int _{\tau _1}^{\tau _2}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n-1}\int _Q \Bigg (\varepsilon \vartheta ^5+\delta \varrho ^2+\delta \frac{1}{\vartheta ^3}\Bigg ) \,\mathrm {d}t\end{aligned}$$

due to (2.13). Finally, due to (2.9)–(2.12),

$$\begin{aligned} \frac{p_M(\varrho , \vartheta )}{\vartheta } +\frac{\varrho e_{M,\delta }(\varrho , \vartheta )}{\vartheta } -\varrho s_{\delta ,M}(\varrho , \vartheta ) \end{aligned}$$

is bounded from below by a negative constant, such that the corresponding term can be bounded from below by $-c(\tau _2-\tau _1)$.

Using (2.1), $\int _Q\varrho \,\mathrm {d}x=M_\varepsilon \le M_0$ and (2.13) the terms (V) and (VIII) can be bounded by

$$\begin{aligned} \int _{\tau _1}^{\tau _2}&\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n-1}+\int _{\tau _1}^{\tau _2}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n-2}\int _Q \varrho |\mathbf{u}|^2\,\mathrm {d}x\,\mathrm {d}t\\&\le \,c\int _{\tau _1}^{\tau _2} \Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n-1}+c(\tau _2-\tau _1)\\ {}&\quad +c\int _{\tau _1}^{\tau _2}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n-1}\int _Q {\overline{\vartheta }}\varrho \log (\vartheta ) \,\mathrm {d}t\\&\le \kappa \int _{\tau _1}^{\tau _2} \Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n}+c_\kappa (\tau _2-\tau _1)\\&\quad +\kappa \int _{\tau _1}^{\tau _2}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n-1}\int _Q \Bigg (\varepsilon \vartheta ^5+\delta \varrho ^2+\delta \frac{1}{\vartheta ^3}\Bigg ) \,\mathrm {d}t, \end{aligned}$$

where $\kappa >0$ is arbitrary. Clearly, (IV) vanishes after taking expectations. On account of (2.9)–(2.12) we have

$$\begin{aligned} \frac{p_M(\varrho , \vartheta )}{\varrho \vartheta }+\frac{e_{M,\delta }(\varrho , \vartheta )}{\vartheta } - s_{M,\delta }(\varrho , \vartheta ) \lesssim 1+\frac{\varrho ^{2/3}}{\vartheta }\le \kappa \varrho ^\Gamma + \kappa \frac{1}{\vartheta ^3}+c_\kappa . \end{aligned}$$

(3.11)

Consequently, the estimate for (VII) is analogous to one for (V) and (VIII) above. We quote from [8, equ. (3.107)]

$$\begin{aligned} \!\!\!\!(II) \le \varepsilon \int _{Q}&\frac{1}{\vartheta ^2}g \Bigg | e_M(\varrho , \vartheta ) + \varrho \frac{\partial e_M (\varrho , \vartheta )}{\partial \varrho } \Bigg | |\nabla \varrho | |\nabla \vartheta | \,\mathrm {d}x\nonumber \\&\le \frac{1}{2} \int _{Q} \Bigg [ \delta \Bigg ( \vartheta ^{\Gamma - 2} + \frac{1}{\vartheta ^3} \Bigg ) |\nabla \vartheta |^2 + \frac{\varepsilon \delta }{\vartheta } \Bigg ( \Gamma \varrho ^{\Gamma - 2} + 2 \Bigg ) |\nabla \varrho |^2 \Bigg ] \ \,\mathrm {d}x \end{aligned}$$

(3.12)

provided we choose $\varepsilon = \varepsilon (\delta ) > 0$ small enough. Finally, we clearly have

$$\begin{aligned} (III)&\le \kappa \int _{\tau _1}^{\tau _2}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n-1}\int _Q \Bigg (\varepsilon \vartheta ^5+\delta \frac{1}{\vartheta ^3}\Bigg ) \,\mathrm {d}t\\ {}&\quad +c_\kappa \int _{\tau _1}^{\tau _2}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n-1}\,\mathrm {d}t\\&\le \kappa \int _{\tau _1}^{\tau _2}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n-1}\int _Q \Bigg (\varepsilon \vartheta ^5+\delta \frac{1}{\vartheta ^3}\Bigg )\,\mathrm {d}x\,\mathrm {d}t\\ {}&\quad +\kappa \int _{\tau _1}^{\tau _2}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n}\,\mathrm {d}t+c_\kappa (\tau _2-\tau _1). \end{aligned}$$

Combing everything and choosing $\kappa$ small enough and noticing that $\delta \frac{1}{\vartheta ^3}\le \sigma _{\varepsilon ,\delta }$ yields

$$\begin{aligned}&{\mathbb {E}}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau _2, \cdot ) \ \,\mathrm {d}x\Bigg ]^n +D {\mathbb {E}}\int _{\tau _1}^{\tau _2}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n-1}\int _{Q}\sigma _{\varepsilon ,\delta }\,\mathrm {d}x\nonumber \\&\quad + D{\mathbb {E}}\int _{\tau _1}^{\tau _2}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n-1}\bigg (\int _{Q} \varepsilon \vartheta ^5\,\mathrm {d}x+\int _{\partial Q}\vartheta ^2\,\mathrm {d}{\mathcal {H}}^2\bigg ) \,\mathrm {d}t\nonumber \\&\le \, {\mathbb {E}}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau _1, \cdot ) \ \,\mathrm {d}x\Bigg ]^n+c(\tau _2-\tau _1) \end{aligned}$$

(3.13)

for all $0\le \tau _1<\tau _2$ with some $D>0$. We obtain in particular

$$\begin{aligned}&{\mathbb {E}}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau _2, \cdot ) \ \,\mathrm {d}x\Bigg ]^n +D \int _{\tau _1}^{\tau _2}{\mathbb {E}}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n}\,\mathrm {d}x\\&\le \, {\mathbb {E}}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau _1, \cdot ) \ \,\mathrm {d}x\Bigg ]^n+C(\tau _2-\tau _1). \end{aligned}$$

Applying the Gronwall lemma from [3, Lemma 3.1] and recalling hypothesis (3.7) we obtain

$$\begin{aligned} \begin{aligned} {\mathbb {E}}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^n&\le \exp (-Dt)\bigg ({\mathbb {E}}\Big [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (0, \cdot ) \ \,\mathrm {d}x\Bigg ]^n-\frac{C}{D}\bigg )+\frac{C}{D}\\ {}&\le E_\infty \end{aligned} \end{aligned}$$

uniformly in $\tau >0$. Using this in (3.13) shows

$$\begin{aligned}&{\mathbb {E}}\int _{0}^{\tau }\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]^{n-1}\bigg (\int _{Q}\sigma _{\varepsilon ,\delta }\,\mathrm {d}x+\int _{Q} \varepsilon \vartheta ^5\,\mathrm {d}x+\int _{\partial Q}\vartheta ^2\,\mathrm {d}{\mathcal {H}}^2\bigg ) \,\mathrm {d}t\\&\quad \le \,E_\infty (1+\tau ) \end{aligned}$$

by possibly enlarging $E_\infty$. $\square$

3.3 Stationary solutions

Based on Proposition 3.2 the method from [17] becomes available and we can construct a stationary solution to (3.1)–(3.4) following the ideas from [3] to which we refer for further details. Different to Section 2.3 we consider stationary solutions sitting on the space of trajectories that are defined on the real line R rather than the interval $[0,\infty )$. We will call them entire stationary solutions. This construction is clearly stronger and hence we obtain also stationary solutions in the sense of Definitions 2.1.

Clearly, Definition 2.1 can be easily modified for solutions $((\Omega ,\mathfrak {F},(\mathfrak {F}_t)_{t\ge -T},{\mathbb {P}}),\varrho ,\vartheta ,\mathbf{u},W)$ being defined on $[-T,\infty )$ for some $T>0$. An entire solution is than an object

$$\begin{aligned} ((\Omega ,\mathfrak {F},(\mathfrak {F}_t)_{t\in R},{\mathbb {P}}),\varrho ,\vartheta ,\mathbf{u},W) \end{aligned}$$

which is a solution on $[-T,\infty )$ for any $T>0$. It takes values in the trajectory space

$$\begin{aligned} {\mathcal {T}}&= {\mathcal {T}}_\varrho \times {\mathcal {T}}_\vartheta \times {\mathcal {T}}_{\mathbf{u}}\times {\mathcal {T}}_W,\\ {\mathcal {T}}_\varrho&=\left( L^2_\mathrm{loc}(R; W^{1,2} (Q; R^d) ), w \right) \cap C_\mathrm{weak, loc}(R; L^\Gamma (Q)),\\ {\mathcal {T}}_{\vartheta }&=\left( L^2_\mathrm{loc}(R; W^{1,2} (Q; R^d) ), w \right) \cap \left( L^\infty _\mathrm{loc}(R; L^{4} (Q) ), w^* \right) \\ {\mathcal {T}}_{\mathbf{u}}&=\left( L^2_\mathrm{loc}(R; W^{1,2}_0 (Q; R^d) ), w \right) ,\quad {\mathcal {T}}_W= C_\mathrm{loc,0}(R; \mathfrak {U}_0 ), \end{aligned}$$

where $C_{\mathrm{loc}, 0}$ denotes the space of continuous functions vanishing at 0. We denote by $\mathfrak {P}({\mathcal {T}})$ the set of Borel probability measures on ${\mathcal {T}}$.

We say that an entire solution to (3.1)–(3.4) of the problem (3.1)–(3.4) is stationary if its law ${\mathcal {L}}_{{\mathcal {T}}}[\varrho ,\vartheta , \mathbf{u}, W]$ is shift invariant in the trajectory space ${\mathcal {T}}$, meaning

$$\begin{aligned} {\mathcal {L}}_{{\mathcal {T}}} \left[ {\mathcal {S}}_\tau [\varrho ,\vartheta , \mathbf{u}, W] \right] = {\mathcal {L}}_{{\mathcal {T}}}[\varrho ,\vartheta , \mathbf{u}, W] \ \text{ for } \text{ any }\ \tau \in R, \end{aligned}$$

with the time shift operator

$$\begin{aligned} {\mathcal {S}}_\tau [\varrho ,\vartheta , \mathbf{u}, W] (t) = [\varrho (t + \tau ),\vartheta (t+\tau ), \mathbf{u}(t + \tau ), W (t + \tau )-W(\tau )], \ t \in R, \ \tau \in R. \end{aligned}$$

Proposition 3.3

Let the assumptions of Theorem 2.1 be valid and let $\varepsilon \le \varepsilon _0$ and $\delta ,{\overline{\vartheta }}>0$ be given. Let

$$\begin{aligned} \big ((\Omega ,\mathfrak {F},(\mathfrak {F}_t)_{t\ge 0},{\mathbb {P}}),\varrho ,\vartheta ,\mathbf{u},W) \end{aligned}$$

be a weak martingale solution of the problem (3.1)–(3.4) (in the sense of Definition 2.1 with the obvious modifications) such that

$$\begin{aligned} \mathrm{ess} \limsup _{t \rightarrow 0 +} {\mathbb {E}} \left[ E_H^{\delta ,{\overline{\vartheta }}}(t)^4 \right] < \infty . \end{aligned}$$

(3.14)

Then there is a sequence $T_n \rightarrow \infty$ and an entire stationary solution

$$\begin{aligned} \big (({\tilde{\Omega }},{\tilde{\mathfrak {F}}},({\tilde{\mathfrak {F}}}_t)_{t\in R},{\tilde{{\mathbb {P}}}}),{{\tilde{\varrho }}},\vartheta ,{\tilde{\mathbf{u}}}, {\tilde{W}}) \end{aligned}$$

such that

$$\begin{aligned} \frac{1}{T_n} \int _0^{T_n} {\mathcal {L}}_{{\mathcal {T}}}\left[ {\mathcal {S}}_t \left[ \varrho ,\vartheta , {\mathbf{u}}, {W} \right] \right] \,\mathrm {d}t\rightarrow {\mathcal {L}}_{{\mathcal {T}}} \left[ {{\tilde{\varrho }}}, {\tilde{\mathbf{u}}},{\tilde{\vartheta }}, {\tilde{W}} \right] \ \text{ narrowly } \text{ as }\ n \rightarrow \infty . \end{aligned}$$

Proof

Let $[\varrho , \mathbf{u}, W]$ be a dissipative martingale solution on $[0,\infty )$ defined on some stochastic basis $(\Omega ,\mathfrak {F},(\mathfrak {F}_{t})_{t\ge 0},{\mathbb {P}})$ and satisfying (3.14). We define the probability measures

$$\begin{aligned} \nu _S \equiv \frac{1}{S} \int _0^S {\mathcal {L}}_{{\mathcal {T}}} (S_t [\varrho ,\vartheta , \mathbf{u}, W]) \,\mathrm {d}t\in \mathfrak {P} ({\mathcal {T}}) . \end{aligned}$$

(3.15)

We tacitly regard functions defined on time intervals $[-t,\infty )$ as trajectories on R by extending them to $s\le -t$ by the value at $-t$. As in [3, Prop. 5.1] we can show that the family of measures $\{\nu _{S};\,S>0\}$ is tight on ${\mathcal {T}}$. In fact, Proposition 3.2 yields

$$\begin{aligned}&{\mathbb {E}} \left[ \sup _{s \in [-T,T]} {\mathcal {E}}_\delta ^m (s+t) \right] + {\mathbb {E}} \left[ \int _{-T\vee -t}^T \int _{Q} \bigg (|\nabla \varrho |^2+|\nabla \vartheta |^2+|\nabla \mathbf{u}|^{2}\bigg )(s+t) \ \,\mathrm {d}x \,\mathrm{d}s \right] \\&\lesssim {\mathbb {E}} \left[ {\mathcal {E}}^{m}(0) \right] +c. \end{aligned}$$

This gives the same bounds on $\varrho$ and $\mathbf{u}$ as in [3] and we control additionally

$$\begin{aligned} {\mathbb {E}} \left[ \sup _{s \in [-T,T]} \bigg (\int _Q\vartheta ^4\,\mathrm {d}x\bigg )^m (s+t) \right] + {\mathbb {E}} \left[ \int _{-T\vee -t}^T \int _{Q} |\nabla \vartheta |^2(s+t) \ \,\mathrm {d}x \right] \end{aligned}$$

which implies tightness of $\frac{1}{S} \int _0^S {\mathcal {L}}_{{\mathcal {T}}} (S_t [\vartheta ]) \,\mathrm {d}t$. Note also that we have control of $\nabla \varrho$ due to $\varepsilon >0$ which is different from [3].

Due to [3, Lemma 5.2], if the narrow limit of

$$\begin{aligned} \nu _{\tau , S_{n}} \equiv \frac{1}{S_{n}} \int _0^{S_{n}} {\mathcal {L}} (S_{t + \tau } [\varrho , \mathbf{u}, W]) \,\mathrm {d}t\end{aligned}$$

in $\mathfrak {P} ({\mathcal {T}})$ as $n \rightarrow \infty$ exists for some $\tau =\tau _0 \in R$ then it exists for all $\tau \in R$ and is independent of the choice of $\tau$. Applying Jakubowski–Skorokhod’s theorem [18], we infer the existence of a sequence $S_{n}\rightarrow \infty$ and $\nu \in \mathfrak {P} ({\mathcal {T}})$ so that $\nu _{0, S_n} \rightarrow \nu$ narrowly in $\mathfrak {P} ({\mathcal {T}})$ as well as $\nu _{\tau , S_n} \rightarrow \nu$ narrowly for all $\tau \in R$. Accordingly, the limit measure $\nu$ is shift invariant in the sense that for every $G\in BC({\mathcal {T}})$ and every $\tau \in R$ we have $\nu (G\circ S_{\tau })=\nu (G).$ To conclude the proof of Theorem 3.3, it remains to show that $\nu$ is a law of an entire solution to (3.1)–(3.4).

First of all, we can argue as in [3, Prop. 5.3] to show that for any $S>0$

$$\begin{aligned} \nu _{S}\equiv \frac{1}{S}\int _{0}^{S}{\mathcal {L}}(S_{t}[\varrho , \mathbf{u}, W])\,\mathrm {d}t\in \mathfrak {P}({\mathcal {T}}) \end{aligned}$$

is a dissipative martingale solution on $(- T, \infty )$, provided $(\varrho ,\vartheta , \mathbf{u}, W)$ is a dissipative martingale solution on $(-T,\infty )$ defined on some probability space $(\Omega ,{\mathcal {F}},{\mathbb {P}})$. The idea is to use that (3.1), (3.2) and (3.4) can be written as measurable mappings on the paths space (see the proof of [3, Prop. 5.3] for how to include the stochastic integral). Unfortunately, this is not true for the quantities hidden in $\sigma _{\varepsilon ,\delta }$ appearing in(3.4) and (3.13). However, they are measurable on a subset, where the laws ${\mathcal {L}}(S_{t}[\varrho , \mathbf{u}, W])$ are supported. Recall from (3.9) that $\sigma _{\varepsilon ,\delta }$ belongs a.s. to $L^1$ in space and locally in time for any solution. This is enough to arrive at the same conclusion.

To finish the proof we argue that the limit $\nu$ is the law of an entire solution to (3.1)–(3.4). Now we consider the measures $\nu _{\tau , S_n - \tau }$, $n=1,2, \dots$, and $\tau >0$. According to the previous considerations, $\nu _{\tau ,S_{n}-\tau }$ is a dissipative martingale solution to (3.1)–(3.4) on $[-\tau ,\infty )$ and the narrow limit as $n\rightarrow \infty$ exists and equals to $\nu$. Now we take a sequence $\tau _m \rightarrow \infty$ and choose a diagonal sequence such that

$$\begin{aligned} \nu _{\tau _m, S_{n (m)} - \tau _m} \rightarrow \nu \quad {\mathrm as} \quad m \rightarrow \infty . \end{aligned}$$

Applying Jakubowski–Skorokhod’s theorem, we obtain a sequence of approximate processes $[{\tilde{\varrho }}_m, \tilde{\mathbf{u}}_m, {\tilde{W}}_m]$ converging a.s. to a process $[{\tilde{\varrho }}, \tilde{\mathbf{u}}, {\tilde{W}}]$ in the topology of ${\mathcal {T}}$. Moreover, the law of $[{\tilde{\varrho }}_m, \tilde{\mathbf{u}}_m, {\tilde{W}}_m]$ is $\nu _{\tau _{m},S_{n(m)}-\tau _{m}}$ and necessarily the law of $[{\tilde{\varrho }}, \tilde{\mathbf{u}}, {\tilde{W}}]$ is $\nu$. By [2, Thm. 2.9.1] it follows that equations (3.1)–(3.4) as well as (3.5) also hold on the new probability space. The limit procedure on this level is quite easy due to the artificial viscosity: By definition of ${\mathcal {T}}_\varrho$ the sequence ${{\tilde{\varrho }}}_n$ is compact on $L^\Gamma$. Moreover, the strong convergence of ${{\tilde{\vartheta }}}_n$ can be proved exactly as in the deterministic existence theory (see [8, Sec. 3.5.3]). This is enough to pass to the limit in all nonlinearities in (3.1), (3.2) and (3.4). The terms in (3.3) and (3.5) which are not compact (those related to the quantity $\sigma _{\varepsilon ,\delta }$) are convex functions and hence can be dealt with by lower-semicontinuity. $\square$

4 Asymptotic limit

In this section we pass to the limit and the artificial viscosity and the artificial pressure respectively. They crucial point is a uniform-in-time estimate, see (4.7) and (4.8) below, which preserves stationarity in the limit. It has to be combined with pressure estimates which differ on both levels. The key ingredient for estimates (4.7) and (4.8) is the non-homogeneous boundary condition for the temperature, cf. (1.5).

4.1 The vanishing viscosity limit

In this section we start with a stationary solution $(\varrho ,\vartheta ,\mathbf{u})$ to (3.1)–(3.4) existence of which is guaranteed by Proposition 3.3. We prove uniform-in-time estimates and pass subsequently to the limits in $\varepsilon$ and $\delta$.

The entropy balance (3.3) yields after taking expectations and using stationarity

$$\begin{aligned} {\mathbb {E}}\int _{Q}&\frac{1}{\vartheta }\Bigg [{\mathbb {S}}(\vartheta ,\nabla \mathbf{u}):\nabla \mathbf{u}+\frac{\kappa _\delta (\vartheta )}{\vartheta }|\nabla \vartheta |^2+\delta \frac{1}{\vartheta ^2}\Bigg ] \,\mathrm {d}x \\&\le {\mathbb {E}}\int _Q\bigg (-2\varepsilon \varrho +2\varepsilon M_\varepsilon \bigg )\frac{1}{\vartheta }\Bigg ( \frac{p_M(\varrho , \vartheta )}{\varrho }+ e_{M,\delta }(\varrho , \vartheta ) - \vartheta s_{M,\delta }(\varrho , \vartheta ) \Bigg )\,\varphi \,\mathrm {d}x\\&\quad +{\mathbb {E}}\int _Q\varepsilon \vartheta ^4 \,\mathrm {d}x+{\mathbb {E}}\int _{\partial Q}\frac{d(\vartheta )}{\vartheta }(\vartheta -\Theta _0)\,\mathrm {d}{\mathcal {H}}^2. \end{aligned}$$

On account of (3.11) the first two terms can be bounded by

$$\begin{aligned} c\,{\mathbb {E}}\Bigg [ \int _{Q} \delta \varrho ^\Gamma \ \,\mathrm {d}x+1\Bigg ]+\frac{1}{4}{\mathbb {E}}\int _{Q}\delta \frac{1}{\vartheta ^3}\,\mathrm {d}x. \end{aligned}$$

The estimate is independent of $\delta$ if we choose $\varepsilon \le \delta$. Similarly, we obtain

$$\begin{aligned} {\mathbb {E}}\int _{\partial Q}\frac{d(\vartheta )}{\vartheta }(\vartheta -\Theta _0)\,\mathrm {d}{\mathcal {H}}^2&\le \,c\, {\mathbb {E}}\int _{\partial Q}\bigg (\vartheta +|\nabla \vartheta |\bigg )+c\\&\le c\,{\mathbb {E}}\Bigg [ \int _{Q} \vartheta ^4 \ \,\mathrm {d}x+1\Bigg ]+\frac{1}{4}{\mathbb {E}}\int _{Q}\frac{\kappa _\delta (\vartheta )}{\vartheta }|\nabla \vartheta |^2 \,\mathrm {d}x \end{aligned}$$

using (1.6) and the trace theorem. In conclusion,

$$\begin{aligned} \!\!{\mathbb {E}}\int _{Q} \frac{1}{\vartheta }\Bigg [{\mathbb {S}}(\vartheta ,\nabla \mathbf{u}):\nabla \mathbf{u}+\frac{\kappa _\delta (\vartheta )}{\vartheta }|\nabla \vartheta |^2+\delta \frac{1}{\vartheta ^2}\Bigg ] \,\mathrm {d}x&\!\lesssim {\mathbb {E}}\Bigg [ \int _{Q} \bigg (\delta \varrho ^\Gamma +\vartheta ^4\bigg ) \ \,\mathrm {d}x\Bigg ]+1\nonumber \\&\!\lesssim {\mathbb {E}}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) \ \,\mathrm {d}x\Bigg ]+1 \end{aligned}$$

(4.1)

independently of $\varepsilon$ and $\delta$ recalling also (2.13) and $\int _Q\varrho \,\mathrm {d}x\le M_\varepsilon \le M_0$. By (2.13) this implies for any $\xi >0$

$$\begin{aligned} {\mathbb {E}}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) \ \,\mathrm {d}x\Bigg ]&\le \,c\, {\mathbb {E}}\Bigg [ \int _{Q} \bigg (\delta \varrho ^\Gamma +\varrho ^{5/3}+\vartheta ^4+\delta \varrho |\log (\vartheta )|\bigg ) \ \,\mathrm {d}x+1\Bigg ]\\&\le \,c {\mathbb {E}}\Bigg [ \int _{Q} \big (\delta \varrho ^\Gamma +\varrho ^{5/3}+\vartheta ^4\big ) \ \,\mathrm {d}x\Bigg ]+\xi \, {\mathbb {E}}\Bigg [ \int _{Q} \delta \frac{1}{\vartheta ^3} \ \,\mathrm {d}x\Bigg ]+c_\xi \\&\le \,c {\mathbb {E}}\Bigg [ \int _{Q} \bigg (\delta \varrho ^\Gamma +\varrho ^{5/3}+\vartheta ^4\bigg ) \ \,\mathrm {d}x\Bigg ]\\&+c\xi \,{\mathbb {E}}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) \ \,\mathrm {d}x\Bigg ]+c_\xi \end{aligned}$$

such that

$$\begin{aligned} {\mathbb {E}}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) \ \,\mathrm {d}x\Bigg ]&\lesssim {\mathbb {E}}\Bigg [ \int _{Q} \bigg (\delta \varrho ^\Gamma +\varrho ^{5/3}+\vartheta ^4\bigg ) \ \,\mathrm {d}x\Bigg ]+1 \end{aligned}$$

(4.2)

independently of $\varepsilon$ and $\delta$.

In (3.5) we choose $n=2$, apply expectations and use stationarity to obtain

$$\begin{aligned}&2\overline{\vartheta } {\mathbb {E}}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]\int _{Q}\sigma _{\varepsilon ,\delta }\,\mathrm {d}x\nonumber \\&\quad + 2{\mathbb {E}}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]\int _{Q} \varepsilon \vartheta ^5\,\mathrm {d}x\nonumber \\&\quad +2{\mathbb {E}}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]\int _{\partial Q}\frac{d}{\vartheta }\Bigg (\vartheta -\overline{\vartheta }\Bigg )(\vartheta -\Theta _0)\,\mathrm {d}{\mathcal {H}}^2\nonumber \\&\quad + 4\varepsilon {\overline{\vartheta }}{\mathbb {E}}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]\int _{Q} \left[ \delta \varrho ^2 + \frac{\delta \Gamma }{\Gamma - 1} \varrho ^{\Gamma } +\frac{1}{2}\varrho |\mathbf{u}|^2\right] \ \,\mathrm {d}x \nonumber \\&\quad +4\varepsilon {\mathbb {E}}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]\int _{Q}\varrho \Bigg ( \frac{p(\varrho , \vartheta )}{\varrho \vartheta } +\frac{ e_\delta (\varrho , \vartheta )}{\vartheta } -s_\delta (\varrho , \vartheta ) \Bigg )\,\mathrm {d}x\nonumber \\&\le 2\varepsilon {\mathbb {E}}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]\int _{Q}\frac{\overline{\vartheta }}{\vartheta ^2}\bigg (e_{M,\delta }(\varrho ,\vartheta )+\varrho \frac{\partial e_M}{\partial \varrho }(\varrho ,\vartheta )\bigg )\nabla \varrho \cdot \nabla \vartheta \,\mathrm {d}x\nonumber \\&\quad + {\mathbb {E}}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]\bigg ( \int _{Q} \sum _{k \ge 1} \varrho | \mathbf{F}_{k,\varepsilon } (\varrho ,\vartheta ,\mathbf{u}) |^2 \ \,\mathrm {d}x \bigg )\nonumber \\&\quad + 2\varepsilon M_\varepsilon {\mathbb {E}}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ] \int _{Q} \left( 2\delta \varrho + \frac{\delta \Gamma }{\Gamma - 1} \varrho ^{\Gamma - 1}+\frac{1}{2}|\mathbf{u}|^2\right) \ \,\mathrm {d}x \nonumber \\&\quad + 2\varepsilon M_\varepsilon {\overline{\vartheta }}{\mathbb {E}}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ] \int _{Q} \left( \frac{p(\varrho , \vartheta )}{\varrho \vartheta }+\frac{e_\delta (\varrho , \vartheta )}{\vartheta } - s_\delta (\varrho , \vartheta ) \right) \ \,\mathrm {d}x \nonumber \\&\quad +\sum _{k=1}^\infty {\mathbb {E}}\bigg ( \int _{Q}\varrho \mathbf{F}_{k,\varepsilon } (\varrho ,\vartheta ,\mathbf{u})\cdot \mathbf{u}\,\mathrm {d}x\bigg )^2\nonumber \\&=:(I)+(II)+(III)+(IV)+(V). \end{aligned}$$

(4.3)

Arguing as in the proof of Proposition 3.2 but paying attention to the $\varepsilon$- and $\delta$-dependence we have

$$\begin{aligned} (I)&\le \frac{1}{4}{\mathbb {E}}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]\int _Q\sigma _{\varepsilon ,\delta }\,\mathrm {d}x,\\ (II)&\le \,c\, {\mathbb {E}}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ],\\ (III)&\le \varepsilon {\overline{\vartheta }}{\mathbb {E}}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]\int _{Q} \left[ \delta \varrho ^2 + \frac{\delta \Gamma }{\Gamma - 1} \varrho ^{\Gamma }\right] \ \,\mathrm {d}x\\&+\,c\, {\mathbb {E}}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]+\frac{1}{4}{\mathbb {E}}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]\int _Q\sigma _{\varepsilon ,\delta }\,\mathrm {d}x,\\ (IV)&\le \varepsilon {\overline{\vartheta }}{\mathbb {E}}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]\int _{Q} \frac{\delta \Gamma }{\Gamma - 1} \varrho ^{\Gamma } \ \,\mathrm {d}x+c\\&\quad +\frac{1}{4}{\mathbb {E}}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]\int _Q\sigma _{\varepsilon ,\delta }\,\mathrm {d}x,\\ (V)&\le \,c\, {\mathbb {E}}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) (\tau , \cdot ) \ \,\mathrm {d}x\Bigg ]. \end{aligned}$$

Again these estimates are also uniform in $\delta$ if we choose $\varepsilon$ small enough compared to $\delta$. Recalling (2.13) and $\int _Q\varrho \,\mathrm {d}x\le M_\varepsilon \le M_0$ we thus obtain

$$\begin{aligned}&\Vert \vartheta \Vert _{L^6(Q)} {\mathop {\sim }\limits ^{<}}\Vert \vartheta \Vert _{W^{1,2}(Q)} {\mathop {\sim }\limits ^{<}}1+ \int _Q \frac{\overline{\vartheta }}{\vartheta } \frac{\kappa (\vartheta )}{\vartheta } |\nabla \vartheta |^2 \mathrm{d}x + \int _{\partial Q}\vartheta ^2 \,\mathrm {d}{\mathcal {H}}^2,\\&{\mathbb {E}} \left[ \int _Q E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ,\mathbf{u}) \ \,\mathrm {d}x\Vert \vartheta \Vert _{L^6(Q)} \right] {\mathop {\sim }\limits ^{>}} {\mathbb {E}} \left[ \Vert \vartheta \Vert _{L^4(Q)}^4 \Vert \vartheta \Vert _{L^6(Q)} \right] - {\mathbb {E}} \left[ \Vert \vartheta \Vert _{L^6(Q)} \right] \\ \qquad \qquad&\ge {\mathbb {E}} \left[ \Vert \vartheta \Vert _{L^{30/7}(Q)}^5 \right] - {\mathbb {E}} \left[ \Vert \vartheta \Vert _{L^6(Q)} \right] . \end{aligned}$$

This, inserted in left-hand side of (4.3) yields

$$\begin{aligned} \begin{aligned} {\mathbb {E}} \left[ \Vert \vartheta \Vert _{L^{30/7}(Q)}^{30/7} \right]&+ {\mathbb {E}} \left[ \int _Q E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ,\mathbf{u}) \ \,\mathrm {d}x\int _Q\sigma _{\varepsilon ,\delta }\,\mathrm {d}x \right] \\&{\mathop {\sim }\limits ^{<}}{\mathbb {E}}\Bigg [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) \ \,\mathrm {d}x\Bigg ]+1. \end{aligned} \end{aligned}$$

(4.4)

independently of $\varepsilon$ and $\delta$ using also (4.1).

In order to close the estimate why apply pressure estimates which now can depend on $\delta$. Let us introduce the so–called Bogovskii operator ${\mathcal {B}}$ enjoying the following properties:

$$\begin{aligned} \begin{aligned} {\mathcal {B}} : L^q_0 (Q)&\equiv \! \left\{ f \in L^q(Q) \ \Bigg | \int _Q f \,\mathrm {d}x\!=\! 0 \right\} \rightarrow W^{1,q}_0(Q, R^d),\ 1< q< \infty ,\\ {{\,\mathrm{div}\,}}{\mathcal {B}}[f]&= f ,\\ \Vert {\mathcal {B}}[f] \Vert _{L^r(Q)}&{\mathop {\sim }\limits ^{<}}\Vert \mathbf{g} \Vert _{L^r(Q; R^d)} \ \text{ if }\ f = {{\,\mathrm{div}\,}}\mathbf{g}, \ \mathbf{g} \cdot \mathbf{n}|_{\partial Q} = 0, \ 1< r < \infty , \end{aligned} \end{aligned}$$

(4.5)

see [14, Chapter 3] or [15]. Arguing as in [4, Section 5] (but replacing $\nabla \Delta ^{-1}$ by the Bogovskii operator ${\mathcal {B}}$) we obtain

$$\begin{aligned}&{\mathbb {E}} \left[ \int _{Q} \left[ p(\varrho ,\vartheta )\varrho + \frac{1}{3}\varrho ^{2} |\mathbf{u}|^2 \right] \ \,\mathrm {d}x \right] \\&=c(M_0) \, {\mathbb {E}} \left[ \int _{Q} \Bigg ( p(\varrho ,\vartheta ) + \frac{1}{3}\varrho |\mathbf{u}|^2 \Bigg ) \ \,\mathrm {d}x \right] \\&\quad - {\mathbb {E}} \left[ \int _{Q} \Bigg ( \varrho \mathbf{u}\otimes \mathbf{u}- \frac{1}{3} \varrho |\mathbf{u}|^2 {\mathbb {I}} \Bigg ) : \nabla {\mathcal {B}}( \varrho -M_\varepsilon ) \ \,\mathrm {d}x \right] \\&\quad + {\mathbb {E}} \left[ \int _{Q} \Big (\frac{4}{3}\mu (\vartheta )+\eta (\vartheta ) \Big ) {{\,\mathrm{div}\,}}\mathbf{u}\ \varrho \ \,\mathrm {d}x\,\mathrm {d}t \right] + {\mathbb {E}} \left[ \int _{Q} \varrho \mathbf{u}\cdot {\mathcal {B}} {{\,\mathrm{div}\,}}(\varrho \mathbf{u}) \ \,\mathrm {d}x \right] \\&\quad + 2 \varepsilon \, {\mathbb {E}} \left[ \int _{Q} \varrho _\varepsilon \mathbf{u}_\varepsilon \cdot {\mathcal {B}}\left[ \varrho _\varepsilon -M_\varepsilon \right] \ \,\mathrm {d}x \right] +\varepsilon \, {\mathbb {E}} \left[ \int _{Q} \varrho _\varepsilon ^2 {{\,\mathrm{div}\,}}\mathbf{u}_\varepsilon \ \,\mathrm {d}x \right] \\&=:(I)+(II)+(III)+(IV)+(V)+(VI). \end{aligned}$$

The terms (II) and (IV)–(VI) can be estimated as in [4] (note that they don’t contain $\vartheta$). In fact, we have by (2.13) and the continuity of $\nabla {\mathcal {B}}$

$$\begin{aligned} (II)&\lesssim {\mathbb {E}}\Vert \sqrt{\varrho }\mathbf{u}\Vert _{L^2_x}\Vert \mathbf{u}\Vert _{L^6_x}\Vert \sqrt{\varrho }\nabla {\mathcal {B}}(\varrho -M_\varepsilon )\Vert _{L^3_x}\\&\lesssim {\mathbb {E}}\int _{Q}E_H(\varrho ,\vartheta ,\mathbf{u})\,\mathrm {d}x\Vert \nabla \mathbf{u}\Vert ^2_{L^2_x}+{\mathbb {E}}\Vert \sqrt{\varrho }\nabla {\mathcal {B}}(\varrho -M_\varepsilon )\Vert _{L^3_x}^2\\&\lesssim \,1+ {\mathbb {E}}\int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ,\mathbf{u}) \ \,\mathrm {d}x \end{aligned}$$

provided $\Gamma$ is large enough. Furthermore, we obtain for some $\alpha \in (0,1)$

$$\begin{aligned} (IV)&\lesssim {\mathbb {E}}\Vert \varrho \mathbf{u}\Vert _{L^2_x}^2\lesssim {\mathbb {E}}\Vert \varrho \Vert _{L^3_x}^2\Vert \mathbf{u}\Vert _{L^6_x}^2\lesssim {\mathbb {E}}\Vert \varrho \Vert _{L^3_x}^2\Vert \nabla \mathbf{u}\Vert _{L^2_x}^2\\&\lesssim {\mathbb {E}}\Vert \varrho \Vert _{L^1_x}^{2\alpha }\Vert \varrho \Vert ^{2(1-\alpha )}_{L^\Gamma }\Vert \nabla \mathbf{u}\Vert _{L^2_x}^2\lesssim {\mathbb {E}}\Vert \varrho \Vert ^{2(1-\alpha )}_{L^\Gamma }\Vert \nabla \mathbf{u}\Vert _{L^2_x}^2\\&\lesssim {\mathbb {E}}\int _{Q}E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ,\mathbf{u})\,\mathrm {d}x\Vert \nabla \mathbf{u}\Vert ^2_{L^2_x}+{\mathbb {E}}\Vert \nabla \mathbf{u}\Vert ^2_{L^2_x}\\&\lesssim \,1+ {\mathbb {E}}\int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ,\mathbf{u}) \ \,\mathrm {d}x \end{aligned}$$

using again (2.13), $\int _Q\varrho \,\mathrm {d}x\le M_\varepsilon \le M_0$, (4.4) as well as (4.1). We can estimate (V) and (VI) along the same lines. Using (2.14) and (2.16) we have for $\Gamma$ large enough

$$\begin{aligned} (III)\lesssim {\mathbb {E}} \left[ \int _{Q} \big (\vartheta ^4+\varrho ^4+|\nabla \mathbf{u}|^2\big ) \ \,\mathrm {d}x \right] \lesssim {\mathbb {E}} \left[ 1+\left( \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ,\mathbf{u}) \ \,\mathrm {d}x \right) \right] \end{aligned}$$

due to (2.13) (4.1). Obviously, the same estimate holds for (I) such that we can conclude

$$\begin{aligned} {\mathbb {E}} \left[ \int _{Q} \big (\varrho ^{\Gamma +1}+\vartheta ^4\varrho +\varrho ^2|\mathbf{u}|^2 \big ) \ \,\mathrm {d}x \right] \lesssim {\mathbb {E}} \left[ 1+ \int _{Q} \big (\varrho ^\Gamma +\vartheta ^4+\varrho |\mathbf{u}|^2\big ) \ \,\mathrm {d}x \right] , \end{aligned}$$

(4.6)

using also (4.2). Obviously, we have

$$\begin{aligned} \int _{Q} \big (\varrho ^\Gamma +\vartheta ^4\big ) \ \,\mathrm {d}x&\le \,\xi \int _{Q} \big (\varrho ^{\Gamma +1}+\vartheta ^{30/7}\big ) \ \,\mathrm {d}x+c(\xi ),\\ {\mathbb {E}} \left[ \int _{Q} \varrho |\mathbf{u}|^2 \ \,\mathrm {d}x \right]&\le {{\tilde{\xi }}} {\mathbb {E}} \left[ \int _{Q} \varrho ^2|\mathbf{u}|^2 \ \,\mathrm {d}x \right] +c({{\tilde{\xi }}}) {\mathbb {E}} \left[ \int _{Q} |\nabla \mathbf{u}|^2 \ \,\mathrm {d}x \right] \\&\le {{\tilde{\xi }}} {\mathbb {E}} \left[ \int _{Q} \varrho ^2|\mathbf{u}|^2 \ \,\mathrm {d}x \right] +c({{\tilde{\xi }}})\int _{Q} \big (\varrho ^\Gamma +\vartheta ^4+1\big ) \ \,\mathrm {d}x \end{aligned}$$

for any $\xi ,{{\tilde{\xi }}}>0$ using again (4.1). Consequently, all terms can be absorbed in the left-hand side and we end up with

$$\begin{aligned} {\mathbb {E}} \left[ \int _{Q} \left[ \varrho ^{\Gamma +1}+\vartheta ^4\varrho +\vartheta ^{30/7} + \varrho ^{2} |\mathbf{u}|^2 \right] \ \,\mathrm {d}x \right] \le \,c \end{aligned}$$

(4.7)

as well as

$$\begin{aligned}&{\mathbb {E}} \left[ \left( \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ,\mathbf{u}) \ \,\mathrm {d}x \right) \int _{Q} \sigma _{\varepsilon ,\delta } \,\mathrm {d}x \right] \le \,c \end{aligned}$$

(4.8)

using (4.1).

Estimates (4.7) and (4.8) are sufficient to pass to the limit in (3.1)–(3.4) arguing as in [1, Section 5] (in fact, one has to combine ideas from [2] and [8]). In the limit $\varepsilon \rightarrow 0$, we obtain the following system.

• Equation of continuity.

  $$\begin{aligned} \begin{aligned} \int _{0}^{\infty }&\int _{Q} \left[ \varrho \partial _t \varphi + \varrho \mathbf{u}\cdot \nabla \varphi \right] \ \,\mathrm {d}x \,\mathrm {d}t= 0 \end{aligned} \end{aligned}$$

  (4.9)

  for any $\varphi \in C^\infty _c((0, \infty ) \times Q)$ ${\mathbb {P}}$-a.s.

• Momentum equation.

  $$\begin{aligned} \begin{aligned} \int _0^\infty \partial _t&\psi \int _{Q} \varrho \mathbf{u}\cdot \varvec{\varphi } \ \,\mathrm {d}x \,\mathrm {d}t+ \int _0^\infty \psi \int _{Q} \varrho \mathbf{u}\otimes \mathbf{u}: \nabla \varvec{\varphi } \ \,\mathrm {d}x \,\mathrm {d}t\\&+ \int _0^\infty \psi \int _{Q} p_\delta (\vartheta ,\varrho ) {{\,\mathrm{div}\,}}\varvec{\varphi } \ \,\mathrm {d}x \,\mathrm {d}t\\&- \int _0^\infty \psi \int _{Q} {\mathbb {S}}(\vartheta ,\nabla \mathbf{u}) : \nabla \varvec{\varphi } \ \,\mathrm {d}x \,\mathrm {d}t=- \int _0^\infty \psi \int _{Q} \varrho {\mathbb {F}}(\varrho ,\vartheta ,\mathbf{u}) \cdot \varvec{\varphi } \ \,\mathrm {d}x \ \mathrm{d}W \end{aligned} \end{aligned}$$

  (4.10)

  for any $\psi \in C^\infty _c((0, \infty ))$, $\varvec{\varphi }\in C^\infty (Q; {\mathbb {R}}^3)$ ${\mathbb {P}}$-a.s.

• Entropy balance.

  $$\begin{aligned} \begin{aligned} - \int _0^\infty&\int _{Q} \left[ \varrho s_\delta (\varrho ,\vartheta ) \partial _t \psi + \varrho s_\delta (\varrho ,\vartheta )\mathbf{u}\cdot \nabla \psi \right] \,\varphi \ \,\mathrm {d}x\mathrm {d}t\\\ge&\int _0^\infty \int _{Q} \frac{1}{\vartheta }\Big [{\mathbb {S}}_\delta (\vartheta ,\nabla \mathbf{u}):\nabla \mathbf{u}+\frac{\kappa _\delta (\vartheta )}{\vartheta }|\nabla \vartheta |^2+\delta \frac{1}{\vartheta ^2}\Big ] \psi \,\mathrm {d}x\,\mathrm {d}t \\&+ \int _0^\infty \int _{Q} \frac{\kappa _\delta (\vartheta )\nabla \vartheta }{\vartheta } \cdot \nabla \psi \ \,\mathrm {d}x\,\mathrm {d}t-\int _0^\infty \psi \int _{\partial Q}\frac{d}{\vartheta }(\vartheta -\Theta _0)\,\mathrm {d}{\mathcal {H}}^2\,\mathrm {d}t\end{aligned} \end{aligned}$$

  (4.11)

  for any $\psi \in C^\infty _c((0, \infty ))$, $\varphi \in C^\infty (Q; {\mathbb {R}}^3)$ ${\mathbb {P}}$-a.s.

• Total energy balance.

  $$\begin{aligned} \begin{aligned} - \int _0^T&\partial _t \psi \bigg (\int _{Q} {{\mathcal {E}}_\delta (\varrho ,\vartheta , \mathbf{u}) } \,\mathrm {d}x\bigg ) \,\mathrm {d}t+\int _0^\infty \psi \int _{\partial Q}d(\vartheta -\Theta _0)\,\mathrm {d}{\mathcal {H}}^2\,\mathrm {d}t\\ =&\,\psi (0) \int _{Q} {{\mathcal {E}}_\delta (\varrho _0,\vartheta _0, \mathbf{u}_0 )} \ \,\mathrm {d}x+\int _0^T\int _{Q}\frac{\delta }{\vartheta ^2}\psi \,\mathrm {d}x\,\mathrm {d}t\\&+ \frac{1}{2} \int _0^T \psi \bigg ( \int _{Q} \sum _{k \ge 1} \varrho | \mathbb {F}_{k} (\varrho ,\vartheta , \mathbf{u}) |^2 \ \,\mathrm {d}x \bigg ) \mathrm{d}t\\&+\int _0^T \psi \int _{Q}\varrho {\mathbf{F}}(\varrho ,\vartheta ,\mathbf{u})\cdot \mathbf{u}\,\mathrm {d}W\,\mathrm {d}x\end{aligned} \end{aligned}$$

  (4.12)

  for any $\psi \in C^\infty _c(0, \infty )$ ${\mathbb {P}}$-a.s.

To summarize, we deduce the following.

Proposition 4.1

Let $\delta >0$ be given. Then there exists a stationary weak martingale solution $[\varrho _\delta ,\vartheta _\delta ,\mathbf{u}_\delta ]$ to (4.9)–(4.12). Moreover, we have the estimates

$$\begin{aligned} {\mathbb {E}} \left[ \Vert \vartheta \Vert _{L^{30/7}(Q)}^{30/7} \right]&+ {\mathbb {E}} \left[ \int _Q E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ,\mathbf{u}) \ \,\mathrm {d}x\int _Q\sigma _{\delta }\,\mathrm {d}x \right] \end{aligned}$$

(4.13)

$$\begin{aligned}&{\mathop {\sim }\limits ^{<}}{\mathbb {E}}\Big [ \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ) \ \,\mathrm {d}x\Big ]+1,\nonumber \\ {\mathbb {E}}\int _Q\sigma _{\delta }\,\mathrm {d}x&\lesssim \,1+ {\mathbb {E}}\int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ,\mathbf{u}) \ \,\mathrm {d}x ,\end{aligned}$$

(4.14)

uniformly in $\delta$, where

$$\begin{aligned} \sigma _{\delta }&=\frac{1}{\vartheta }\Big [{\mathbb {S}}(\vartheta ,\nabla \mathbf{u}):\nabla \mathbf{u}+\frac{\kappa (\vartheta )}{\vartheta }|\nabla \vartheta |^2+\frac{\delta }{2}\Big (\vartheta ^{\Gamma -1}+\frac{1}{\vartheta ^2}\Big )|\nabla \vartheta |^2+\delta \frac{1}{\vartheta ^2}\Big ]. \end{aligned}$$

Corollary 4.1

The solution from Proposition 4.1 satisfies the equation of continuity in the renormalised sense.

4.2 The vanishing artificial pressure limit

Though (4.13) and (4.14) are uniform in $\delta$, the final estimates (4.7) and (4.8) are not. Again we have to close the estimate by some pressure bounds. Let $(\varrho ,\vartheta ,\mathbf{u})$ be a stationary solution to (4.9)–(4.12) as obtained in Proposition 4.1. Arguing as in [4, Section 6] (replacing again $\nabla \Delta ^{-1}$ by the Bogovskii operator ${\mathcal {B}}$) we have

$$\begin{aligned} \begin{aligned}&{\mathbb {E}} \left[ \int _{Q} \bigg [ p_\delta (\varrho ,\vartheta ) \varrho ^{\alpha } + \varrho _{\delta }^{1 + \alpha } |\mathbf{u}|^2 \bigg ] \ \,\mathrm {d}x \right] \\&\le c(M_0) \left( {\mathbb {E}} \left[ \int _{Q} \left[ \frac{1}{2} \varrho |\mathbf{u}|^2 + p_\delta (\varrho ,\vartheta ) \right] \ \,\mathrm {d}x \right] + 1 \right) \\&\qquad + {\mathbb {E}} \left[ \int _{Q} \left( \frac{4}{3}\mu (\vartheta )+\eta (\vartheta ) \right) {{\,\mathrm{div}\,}}\mathbf{u}\ \varrho ^\alpha \ \,\mathrm {d}x \right] \\&\qquad + {\mathbb {E}} \left[ \int _{Q} \left( \varrho \mathbf{u}\otimes \mathbf{u}- \frac{1}{3} \varrho |\mathbf{u}|^2 {\mathbb {I}} \right) : \nabla {\mathcal {B}} \left[ \varrho ^\alpha \right] \ \,\mathrm {d}x \right] \\&\qquad + {\mathbb {E}} \left[ \int _{Q} \varrho \mathbf{u}\cdot {\mathcal {B}} [ {{\,\mathrm{div}\,}}(\varrho ^\alpha \mathbf{u}) + (\alpha - 1) \varrho ^\alpha {{\,\mathrm{div}\,}}\mathbf{u}] \ \,\mathrm {d}x \right] \\&=:(I)+(II)+(III)+(IV), \end{aligned} \end{aligned}$$

(4.15)

where $\alpha >0$ will be chosen sufficiently small. As in the proof of [4, Prop. 6.1] we obtain

$$\begin{aligned}&(I)+(III)+(IV)\lesssim {\mathbb {E}} \left[ \int _Q E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ,\mathbf{u}) \ \,\mathrm {d}x\int _Q|\nabla \mathbf{u}|^2\,\mathrm {d}x \right] \\&\quad + {\mathbb {E}} \left[ \int _Q E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ,\mathbf{u}) \ \,\mathrm {d}x \right] +1. \end{aligned}$$

Also we see that

$$\begin{aligned} (III)&\lesssim {\mathbb {E}} \left[ \int _Q|\nabla \mathbf{u}|^2\,\mathrm {d}x \right] + {\mathbb {E}} \left[ \int _Q\big (\varrho ^\gamma +\vartheta ^4\big )\,\mathrm {d}x \right] +1\\ {}&\lesssim {\mathbb {E}} \left[ \int _Q|\nabla \mathbf{u}|^2\,\mathrm {d}x \right] + {\mathbb {E}} \left[ \int _Q E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ,\mathbf{u}) \ \,\mathrm {d}x \right] +1 \end{aligned}$$

choosing $\alpha$ small enough and using (2.14) and (2.16). Combining these estimate with (4.13) and (4.14) we conclude

$$\begin{aligned} \begin{aligned} {\mathbb {E}}\bigg [ \int _Q \big (\delta \varrho ^{\Gamma +\alpha }+&\varrho ^{\gamma +\alpha }+\vartheta ^4\varrho ^\alpha +\varrho ^{1+\alpha }|\mathbf{u}|^2 \big )\,\mathrm {d}x\bigg ]\\ {}&\lesssim {\mathbb {E}} \left[ 1+ \int _{Q} \big (\delta \varrho ^\Gamma +\varrho ^\gamma +\vartheta ^4+\varrho |\mathbf{u}|^2\big ) \ \,\mathrm {d}x \right] , \end{aligned} \end{aligned}$$

(4.16)

recalling also (4.2). As in the proof of (4.7) and (4.8) we deduce

$$\begin{aligned}&{\mathbb {E}} \left[ \int _{Q} \left[ \varrho ^{\gamma +1}+\vartheta ^4\varrho +\vartheta ^{30/7} + \varrho ^{2} |\mathbf{u}|^2 \right] \ \,\mathrm {d}x \right] \le \,c \end{aligned}$$

(4.17)

$$\begin{aligned}&{\mathbb {E}} \left[ \left( \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ,\mathbf{u}) \ \,\mathrm {d}x \right) \int _{Q} |\nabla \mathbf{u}|^2 \,\mathrm {d}x \right] \le \,c, \end{aligned}$$

(4.18)

$$\begin{aligned}&{\mathbb {E}} \left[ \left( \int _{Q} E_H^{\delta ,\overline{\vartheta }}(\varrho ,\vartheta ,\mathbf{u}) \ \,\mathrm {d}x \right) \int _{Q} \sigma _\delta \,\mathrm {d}x \right] \le \,c, \end{aligned}$$

(4.19)

using (4.14). With estimates (4.17) and (4.18) at hand we can follow the lines of [1, Section 6] to pass to the limit $\delta \rightarrow 0$ in (4.9) and (4.10). The limit in (4.9) and (4.10) can be performed as in [1, Section 7] due to (4.17) and (4.19). This finishes the proof of Theorem 2.1.