On the Use and Misuse of the Oberbeck–Boussinesq Approximation

Abstract

The Oberbeck–Boussinesq approximation is the most commonly employed theoretical scheme for the study of natural or mixed convection flows. However, the misunderstanding of this approximated framework is a possibility that may cause the emergence of paradoxes or, at least, incorrect conclusions. In this paper, the basic features of the Oberbeck–Boussinesq approximation are briefly recalled and three simple examples where this theoretical scheme may be misused are provided. Such misuses of the approximation lead to erroneous conclusions that, in the examples presented in this note, entail violations of the principle of mass conservation. A discussion about the Oberbeck–Boussinesq approximation as an asymptotic theory obtained by letting the product of the thermal expansion coefficient and the reference temperature difference tend to zero is also presented.

1. Introduction

The study of natural or mixed convection flows either in fluids or in fluid-saturated porous media is, with a limited number of exceptions, modelled theoretically by claiming the validity of the Oberbeck–Boussinesq approximation. There are several thorough and comprehensive analyses of how this approximation can be established starting from a general formulation of the local balance equations of mass, momentum and energy for a fluid. Beyond the many textbooks of fluid dynamics and convection heat transfer, we mention the analyses of this topic presented in Chapter 8 of the book by Zeytounian [1], as well as in the papers by Rajagopal, Ruzicka, and Srinivasa [2] and by Zeytounian [3]. Such discussions on the origin and on the range of applicability of the Oberbeck–Boussinesq approximation stem from pioneering papers such as those by Spiegel and Veronis [4], by Gray and Giorgini [5] and by Hills and Roberts [6]. On the other hand, further papers offer interesting perspectives on possible extended forms of the Oberbeck–Boussinesq approximation [7,8,9,10,11,12].

The aim of this short paper is neither a comprehensive analysis of the Oberbeck–Boussinesq approximation nor a novel deduction of the approximate mathematical model. For the former purpose, we refer to the extremely large body of literature on the topic, where the few references cited above are just a definitely limited sample, while the books by Gebhart, Jaluria, Mahajan, and Sammakia [[13] and by Straughan [14] provide comprehensive surveys. For the latter purpose, we mention the recent paper by Barletta [15]. In fact, the objective of this paper is to highlight an aspect of the Oberbeck–Boussinesq approximation that may be the source of pitfalls, i.e., the duality of variable fluid density and constant fluid density. If it is recognised that the approximation scheme predicates a variable density which is pressure independent and varying with the temperature, in some instances one may forget that such variable density serves only to define the buoyancy force within the local momentum balance equation. Utilising the variable density outside this very specific context may lead to unphysical predictions and, hence, to incorrect conclusions. Such conclusions are incorrect as they usually conflict with the principle of mass conservation. A final discussion about the interpretation of the Oberbeck–Boussinesq approximation as a limiting case of the general local balance equations of mass, momentum and energy is presented.

2. A Minimalistic Survey of the Oberbeck–Boussinesq Model

There are several detailed and thorough descriptions of the Oberbeck–Boussinesq model for natural and mixed convection flows. Such a model describes flow conditions with a weakly variable density due to temperature changes in the fluid.

2.1. The Buoyant Flow Governing Equations

The basic system of partial differential equations, expressing the local mass, momentum and energy balances, is given by

(1)$\begin{array}{c}\mathbf{\nabla }\mathbf{·}\mathbf{u}=0,\end{array}$

(2)$\begin{array}{c}{\rho }_0[\frac{\partial \mathbf{u}}{\partial t}+(\mathbf{u}\mathbf{·}\mathbf{\nabla })\mathbf{u}]=-\mathbf{\nabla }P+\rho \mathbf{g}+\mu {\nabla }^2\mathbf{u},\end{array}$

(3)$\begin{array}{c}\frac{\partial T}{\partial t}+(\mathbf{u}\mathbf{·}\mathbf{\nabla })T=\alpha {\nabla }^{2}T,\end{array}$

where u is the velocity, P is the pressure, T is the temperature, t is time, $\mu$ is the dynamic viscosity, g is the gravitational acceleration, and $\alpha$ is the thermal diffusivity. There are two densities in Equations (1)–(3). One is the reference density, ${\rho }_0$, i.e., the fluid density evaluated at the reference temperature, $T_0$. On the other hand, $\rho$ denotes the fluid density evaluated at the local temperature T through the linear equation of state,

(4)$\begin{array}{c}\rho ={\rho }_0[1-\beta (T-T_0)],\end{array}$

where $\beta$ is the thermal expansion coefficient. We mention that the linear equation of state (4) is to be replaced by a quadratic equation of state in special situations such as pure water close to 4 °C.

It must be stressed that, in Equations (1)–(4), ${\rho }_0, \mu , \alpha ,$ and $\beta$ are constant, pressure independent, fluid properties evaluated at $T_0$. In order to maximise the reliability of the approximation, a judicious choice of the reference temperature is an average temperature over the flow spatial domain and over the time interval of the flow process. We also implicitly assume, for the local energy balance (3), that the viscous heating effect is negligible and that no internal heat source is present.

2.2. The Buoyant Seepage Flow in a Porous Medium

When the Oberbeck–Boussinesq scheme is to be applied to the seepage flow in a porous medium, then the local mass balance Equation (1) is rewritten as

(5)$\begin{array}{c}\mathbf{\nabla }\mathbf{·}\mathbf{U}=0.\end{array}$

Here, U denotes the seepage, or Darcy’s, velocity in the porous medium. The seepage velocity results from a representative-elementary-volume average of the pore-scale fluid velocity, u. Furthermore, the local momentum balance equation may be modelled through Darcy’s law, namely,

(6)$\begin{array}{c}\frac{\mu }{K}\mathbf{U}=-\mathbf{\nabla }P+\rho \mathbf{g},\end{array}$

where K is the permeability of the medium. Then, Equation (6) supersedes Equation (2), in this case. Also Equation (3) is to be reformulated when one studies the seepage flow in a porous medium. In fact, the local energy balance reads

(7)$\begin{array}{c}\sigma \frac{\partial T}{\partial t}+(\mathbf{U}\mathbf{·}\mathbf{\nabla })T={\alpha }_m{\nabla }^{2}T,\end{array}$

where $\sigma$ is the ratio between the average heat capacity of the saturated porous medium and the heat capacity of the fluid, while ${\alpha }_m$ is the average thermal diffusivity of the saturated porous medium.

2.3. Hydrostatic Pressure, Buoyancy Force

The gravitational force term, $\rho \mathbf{g}$ in either Equation (2), or in Equation (6), is the one and only quantity where the use of Equation (4) is allowed. Thus, we can write

(8)$\begin{array}{c}-\mathbf{\nabla }P+\rho \mathbf{g}=-\mathbf{\nabla }\left[P-\underset{{\displaystyle P_h}}{\underbrace{(P_e+{\rho }_0\mathbf{g}\mathbf{·}\mathbf{x})}}\right]+\mathbf{b},\\ -\mathbf{\nabla }P+\rho \mathbf{g}=-\mathbf{\nabla }(P-P_h)+\mathbf{b},\end{array}$

where

(9)$\begin{array}{c}{P}_h=P_e+{\rho }_0 \mathbf{g}\mathbf{·}\mathbf{x}\end{array}$

is the hydrostatic pressure, with $\mathbf{x}=(x,y,z)$ denoting the position vector, and

(10)$\begin{array}{c}\mathbf{b}=-{\rho }_0\beta (T-T_0)\mathbf{g}\end{array}$

is the buoyancy force. Both P and $P_h$ depend on the spatial position. It is a common practice to use a specific notation for their local difference,

(11)$\begin{array}{c}p=P-P_h.\end{array}$

We call p hydrodynamic pressure, for the sake of brevity, to mark the distinction from the pressure and the hydrostatic pressure. In Equation (9), $P_e$ is an arbitrary constant value which, physically, can be assigned in order to fix a possible overall external pressurisation of the fluid. In most usual cases, $P_e$ is assumed to coincide with the atmospheric pressure.

Out of the present context, where the hydrostatic pressure and the buoyancy force are defined on the basis of a chosen reference temperature, $T_0$, the fluid density is considered to be a constant equal to ${\rho }_0$. If one violates this simple rule, then erroneous conclusions could be drawn which will typically lead to violations of the mass conservation principle.

3. The Rectangular Cavity with Side Heating

Let us consider the classical problem of two–dimensional natural convection in a rectangular cavity with rigid and impermeable side boundaries kept at uniform, but different, temperatures, $T_1$ and $T_2$, while the upper and lower sides are kept adiabatic. It is not restrictive to assume $T_1>T_2$. A sketch of the system is provided in Figure 1.

It is well–known that, with sufficiently small differences $T_1$−$T_2$, a steady-state natural convection flow occurs in the cavity. In non-dimensional terms, this restriction is equivalent to a sufficiently small Rayleigh number. The steady-state flow is cellular in character, with one or more convective cells. For a given fluid, the number of stationary convective cells depends on the aspect ratio, $H/L$, where H and L are the half–height and half–width of the cavity, and on the Rayleigh number.

We denote with S the horizontal midplane surface with $y=0$ and $x\in [-L,L]$. Since the flow system is modelled as two-dimensional in the xy plane, S is considered with unit depth in the z direction. Then, we can immediately conclude that

(12)$\begin{array}{c}\underset{S}{\int }v\mathrm{d}S=0\hspace{1em}\mathrm{or}\hspace{1em}\underset{-L}{\overset{L}{\int }}v\mathrm{d}x=0,\end{array}$

as a straightforward consequence of Equation (1). In fact, one can just integrate Equation (1) over the upper half-domain $\{x\in [-L,L],y\in [0,H\left]\right\}$ and then employ Gauss’ theorem by recalling that S is the only permeable boundary of such a domain. Here, v is the y component of u. Then, we may wonder how we could evaluate the mass flow rate, $\dot{m}$, across S. There is a correct way and an incorrect way. The correct way is by applying the principle that the fluid density is to be considered constant and equal to ${\rho }_0$, so that Equation (12) yields

(13)$\begin{array}{c}\dot{m}={\rho }_0\underset{-L}{\overset{L}{\int }}v\mathrm{d}x=0.\end{array}$

The incorrect way is by employing the variable density given by Equation (4), so that Equation (12) yields

(14)$\begin{array}{c}\dot{m}=\underset{-L}{\overset{L}{\int }}\rho v\mathrm{d}x={\rho }_0\underset{-L}{\overset{L}{\int }}v\mathrm{d}x-{\rho }_0\beta \underset{-L}{\overset{L}{\int }}(T-T_0)v\mathrm{d}x=-{\rho }_0\beta \underset{-L}{\overset{L}{\int }}(T-T_0)v\mathrm{d}x<0.\end{array}$

The difference between the correct method expressed by Equation (13) and the incorrect method expressed by Equation (14) is that, in the former, we employ the reference density, ${\rho }_0$, while, in the latter, we employ the variable density, $\rho$. The reason why Equation (14) yields $\dot{m}<0$ is that, at least with a sufficiently small temperature difference, $T_1-T_2$, both v and $T-T_0$ are odd functions of x with $x\in [-L,L]$, positive for $x\in (-L,0)$ and negative for $x\in (0,L)$ (see, for instance, de Vahl Davis [16]). Thus, their product, $(T-T_0)v$, is an even and positive function of x throughout the domain of integration $x\in [-L,L]$, so that Equation (14) leads to the conclusion that $\dot{m}<0$. The numerical solution discussed by de Vahl Davis [16] is complemented by a straightforward analytical solution, which predicts the same symmetry for v and $T-T_0$ and which holds for a very tall cavity, $H\gg L$. Such an asymptotic solution is briefly outlined in Appendix A with an explicit, though incorrect, evaluation of the mass flow rate based on the formula given by Equation (14).

We mention that Equation (12) holds also for any other $y=\mathrm{constant}$ plane S, as a consequence of Equation (1). We also note that, here, our arguments are strictly relative to the example examined in this section. Certainly, there are situations where the mass flow rate, evaluated coherently with the Oberbeck–Boussinesq approximation, is nonzero. For instance, this is the case for the flow system to be discussed in Section 4.

We said that Equation (14) expresses the incorrect way to evaluate $\dot{m}$, since the conclusion $\dot{m}\ne 0$ is an evident violation of the mass conservation within the upper half–domain $\{x\in [-L,L],y\in [0,H\left]\right\}$ or in the lower half–domain, $\{x\in [-L,L],y\in [-H,0\left]\right\}$. In fact, both for the upper and the lower half–domains, S would be the only permeable boundary and it would be crossed by a net mass flow rate. Such a situation, in a stationary regime, yields a violation of the principle of mass conservation.

4. Mixed Convection Duct Flow

Let us consider the internal mixed convection in a duct with an impermeable wall having an increasing temperature along the streamwise direction. Such a behaviour is observed, for instance, when the duct wall is subject to an incoming uniform heat flux. We consider the case of stationary flow.

As sketched in Figure 2, we consider the region V delimited by the cross–sections, $S_1$ and $S_2$. Since, the fluid is heated in the streamwise direction, we have an average temperature evaluated at the cross–section $S_1$, denoted by $T_1$, smaller than the average temperature evaluated at the cross–section $S_2$, denoted by $T_2$. Then, a judicious choice of the reference temperature, $T_0$, for the Oberbeck–Boussinesq approximation in the domain V is the volume-averaged temperature,

(15)$\begin{array}{c}{T}_0=\underset{V}{\int }T\mathrm{d}V.\end{array}$

Such a volume–averaged temperature value, $T_0$, is larger than $T_1$ and smaller than $T_2$. Let us evaluate the average velocities across $S_1$ and $S_2$,

(16)$\begin{array}{c}{u}_{m1}=\frac{1}{S_1}\underset{S_1}{\int }\mathbf{u}\mathbf{·}{\hat{\mathbf{n}}}_1\mathrm{d}S,\hspace{1em}{u}_{m2}=\frac{1}{S_2}\underset{S_2}{\int }\mathbf{u}\mathbf{·}{\hat{\mathbf{n}}}_2\mathrm{d}S,\end{array}$

where ${\hat{\mathbf{n}}}_1$ and ${\hat{\mathbf{n}}}_2$ are the unit vectors of the surfaces $S_1$ and $S_2$, oriented in the streamwise direction as shown in Figure 2. Equation (1), after an integration over V and use of Gauss’ theorem, yields the equality $S_1{u}_{m1}=S_2{u}_{m2}$. One can wonder how one can evaluate the mass flow rates across $S_1$ and $S_2$, i.e., ${\dot{m}}_1$ and ${\dot{m}}_2$, respectively. The right way is by assuming that the fluid density is to be considered constant and equal to ${\rho }_0$ all over V, so that

(17)$\begin{array}{c}{\dot{m}}_1={\rho }_0\underset{S_1}{\int }\mathbf{u}\mathbf{·}{\hat{\mathbf{n}}}_1\mathrm{d}S,\hspace{1em}{\dot{m}}_2={\rho }_0\underset{S_2}{\int }\mathbf{u}\mathbf{·}{\hat{\mathbf{n}}}_2\mathrm{d}S.\end{array}$

Since $S_1{u}_{m1}=S_2{u}_{m2}$, Equations (16) and (17) allow one to conclude that ${\dot{m}}_1={\dot{m}}_2$, which is in perfect agreement with the principle of mass conservation.

The incorrect way to evaluate ${\dot{m}}_1$ and ${\dot{m}}_2$ is by employing the variable density, $\rho$, defined by Equation (4). In this case, one would write

(18)$\begin{array}{c}{\dot{m}}_1={\rho }_0\underset{S_1}{\int }\mathbf{u}\mathbf{·}{\hat{\mathbf{n}}}_1\mathrm{d}S-{\rho }_0\beta \underset{S_1}{\int }(T-T_0)\mathbf{u}·{\hat{\mathbf{n}}}_1\mathrm{d}S,\\ {\dot{m}}_2={\rho }_0\underset{S_2}{\int }\mathbf{u}·{\hat{\mathbf{n}}}_2\mathrm{d}S-{\rho }_0\beta \underset{S_2}{\int }(T-T_0)\mathbf{u}·{\hat{\mathbf{n}}}_2\mathrm{d}S.\end{array}$

From Equation (18) and from the equality $S_1{u}_{m1}=S_2{u}_{m2}$, one has

(19)$\begin{array}{c}{\dot{m}}_1-{\dot{m}}_2=-{\rho }_0\beta (\underset{S_1}{\int }T\mathbf{u}\mathbf{·}{\hat{\mathbf{n}}}_1\mathrm{d}S-\underset{S_2}{\int }T\mathbf{u}\mathbf{·}{\hat{\mathbf{n}}}_2\mathrm{d}S).\end{array}$

By employing again the equality $S_1{u}_{m1}=S_2{u}_{m2}$ and by recalling that the average temperature over $S_1$ is smaller than the average temperature over $S_2$, one can immediately conclude that ${\dot{m}}_1-{\dot{m}}_2>0$, i.e., the mass rate flowing through $S_1$ is larger than the mass rate flowing through $S_2$ which, in a stationary regime, means a violation of the principle of mass conservation.

5. A Vertical Porous Slab Separating Two Fluid Reservoirs

Let us consider an infinitely wide wall separating two fluid reservoirs kept at different temperatures, $T_1$ and $T_2$. As shown in Figure 3, the wall has a porous insertion bounded by two planes $S_1$ and $S_2$. The fluid in the left-hand reservoir is the same as that in the right-hand reservoir and they are both in a rest state. The same fluid saturates the porous slab. The Oberbeck–Boussinesq approximation can be applied by assuming the reference temperature, $T_0=(T_1+T_2)/2$.

One can evaluate the pressure distribution on the boundaries $S_1$ and $S_2$ as the hydrostatic pressure. By employing Equation (9), the pressure distribution on the plane $S_1$ is given by

(20)$\begin{array}{c}{P}_1=P_{e1}-{\rho }_0 g y,\end{array}$

while, on $S_2$, we have

(21)$\begin{array}{c}{P}_2=P_{e2}-{\rho }_0 g y.\end{array}$

Equations (20) and (21) have been used to formulate the pressure conditions at the porous slab boundaries in Refs. [17,18,19,20]. The pressurisation constants, $P_{e1}$ and $P_{e2}$, can be equal or not depending on the conditions externally imposed on the two reservoirs. For instance, one may have a situation where the reservoir at temperature $T_1$ is compressed on its top, while the reservoir at temperature $T_2$ is open to the atmospheric pressure. In such a situation, one has $P_{e1}>P_{e2}$. Equations (20) and (21) reveal that the pressure distribution on $S_1$ and $S_2$ is y–dependent, but the pressure difference across the porous slab is a constant,

(22)$\begin{array}{c}∆P=P_1-P_2=P_{e1}-P_{e2}.\end{array}$

If none of the two reservoirs is pressurised, then $∆P=0$. If $∆P\ne 0$, one may have a seepage throughflow across the porous slab [18].

As sketched in Figure 3, the porous insertion is bounded above and below by an impermeable material, so that application of Equation (1) yields

(23)$\begin{array}{c}\underset{S_1}{\int }U{d}S=\underset{S_2}{\int }U\mathrm{d}S,\end{array}$

where U is the x component of the seepage velocity, U, across the porous material. By recalling that the fluid density is to be intended as uniform with the value ${\rho }_0$, Equation (23) implies that the mass flow rates across $S_1$ and $S_2$ are equal as required by mass conservation,

(24)$\begin{array}{c}\dot{m}={\rho }_0\underset{S_1}{\int }U\mathrm{d}S={\rho }_0\underset{S_2}{\int }U\mathrm{d}S.\end{array}$

This is the correct and consistent application of the Oberbeck–Boussinesq scheme to this sample case.

There is always the possibility of introducing errors, as in the previous examples. One may enforce the validity of Equation (4) in the evaluation of the hydrostatic pressures in the two reservoirs so that Equations (20) and (21) are replaced by

(25)$\begin{array}{c}{P}_1=P_{e1}-{\rho }_0 g y (1-\frac{\delta }{2}),\end{array}$

and

(26)$\begin{array}{c}{P}_2=P_{e2}-{\rho }_0 g y (1+\frac{\delta }{2}),\end{array}$

where

(27)$\begin{array}{c}\delta =\beta (T_1-T_2).\end{array}$

Equations (25) and (26) yield an y–dependent pressure difference across the porous layer,

(28)$\begin{array}{c}\mathrm{∆}P=P_1-P_2=P_{e1}-P_{e2}+{\rho }_0 \delta g y.\end{array}$

Even in the absence of pressurisation in one of the two reservoirs $(P_{e1}=P_{e2})$, Equation (28) entails a pressure difference across the porous layer which, in turn, leads to the prediction of a horizontal throughflow across the porous layer [21]. Such a stationary throughflow is unphysical as it would lead to a violation of mass conservation. In fact, consistently with the assumption of a variable density, as given by Equation (4), the mass flow rate across $S_1$ is expressed as

(29)$\begin{array}{c}{\dot{m}}_1={\rho }_0(1-\frac{\delta }{2})\underset{S_1}{\int }U{d}S,\end{array}$

while the mass flow rate across $S_2$ is given by

(30)$\begin{array}{c}{\dot{m}}_2={\rho }_0(1+\frac{\delta }{2})\underset{S_2}{\int }U{d}S.\end{array}$

The violation of mass conservation is quantified by a relative error,

(31)$\begin{array}{c}\frac{{\dot{m}}_2-{\dot{m}}_1}{({\dot{m}}_2+{\dot{m}}_1)/2}=\delta ,\end{array}$

where Equations (23), (29) and (30) have been employed. Interestingly enough, should $\delta$ be considered as negligible, $\delta \ll 1$, then the pressure distributions given by Equations (25) and (26), and employed in the study carried out by Vynnycky and Mitchell [21], would turn out to match perfectly the correct pressure distributions expressed by Equations (20) and (21). More precisely, the corrective terms introduced in Equations (25) and (26), i.e., $\mp \delta /2$, can be taken into account consistently only if the solenoidal constraint for the velocity, $\mathbf{\nabla }\mathbf{·}\mathbf{U}=0$, is replaced by the variable-density local mass balance equation, $\mathbf{\nabla }\mathbf{·}\left(\rho \mathbf{U}\right)=0$. An intermediate approximation, such as that used in Ref. [21], where Equations (25) and (26) are employed in combination with $\mathbf{\nabla }\mathbf{·}\mathbf{U}=0$, would be flawed as it leads to a violation of the principle of mass conservation.

6. The Froude Number and the Rayleigh Number

Let us reconsider Equations (8)–(10). The pressure gradient and gravitational force contributions to the local momentum balance are given by the force per unit volume,

(32)$\begin{array}{c}\mathbf{F}=-\mathbf{\nabla }P-{\rho }_0 g {\hat{\mathbf{e}}}_y+{\rho }_0\beta g(T-T_0){\hat{\mathbf{e}}}_y,\end{array}$

where we are assuming that the y–axis is vertical and orientated upward. The unit vector along the y–axis is denoted by ${\hat{\mathbf{e}}}_y$. We can rescale F in order to obtain a dimensionless formulation of the local momentum balance equation. The dimensionless F denoted with an asterisk is given by

(33)$\begin{array}{c}{\mathbf{F}}^{*}=\frac{\mathbf{F}}{\mathrm{∆}{P}_r/L},\end{array}$

where the constant $\mathrm{∆}{P}_r$ is a reference pressure difference and the constant L is a reference length. One can also define the dimensionless coordinates pressure and temperature as

(34)$\begin{array}{c}(x^{*},y^{*},z^{*})=\frac{(x,y,z)}{L},\hspace{1em}{P}^{*}=\frac{P}{{\rho }_0 g L},\hspace{1em}{T}^{*}=\frac{T-T_0}{\mathrm{∆}{T}_r},\end{array}$

where $\mathrm{∆}{T}_r$ is a reference temperature difference. On account of Equations (33) and (34), Equation (32) yields

(35)$\begin{array}{c}{\mathbf{F}}^{*}=-\frac{1}{\mathrm{Fr}}({\mathbf{\nabla }}^{*}{P}^{*}+{\hat{\mathbf{e}}}_y)+\frac{\delta }{\mathrm{Fr}}{T}^{*}{\hat{\mathbf{e}}}_y,\end{array}$

where the Froude number is defined as

(36)$\begin{array}{c}\mathrm{Fr}=\frac{\mathrm{∆}{P}_r}{{\rho }_{0}gL},\end{array}$

and $\delta$ is given by

(37)$\begin{array}{c}\delta =\beta \mathrm{∆}{T}_r.\end{array}$

It is noteworthy that the usual definition is $\mathrm{Fr}=U_0^2/\left(gL\right)$ (see, for instance, the paper by Mayeli and Sheard [10]), where $U_0$ is a reference velocity. In fact, Equation (36) matches perfectly the usual definition provided that one chooses $U_0=\sqrt{\mathrm{∆}{P}_r/{\rho }_0}$. One may also mention that, according to other authors (e.g., Zeytounian [1]), $\mathrm{Fr}$ should be defined as a velocity ratio, so that it is given by the square root of the quantity identified by Mayeli and Sheard [10] as the Froude number. In our discussion, we will rely on the definition given by Equation (36).

It is quite evident that one can, legitimately, set the reference pressure difference $\mathrm{∆}{P}_r={\rho }_{0}gL$. This choice, from Equation (36), would yield $\mathrm{Fr}=1$. In this case, assuming $\delta \to 0$ would justify an assumption of constant density in the mass balance equation, but it would also unavoidably imply a negligible buoyancy contribution to ${\mathbf{F}}^{*}$ as given by Equation (35). Indeed, in this case, the limit $\delta \to 0$ yields

(38)$\begin{array}{c}{\mathbf{F}}^{*}=-{\mathbf{\nabla }}^{*}{P}^{*}-{\hat{\mathbf{e}}}_y.\end{array}$

In Equation (38), there is no buoyancy force, so that the dynamics predicted is that of a perfectly incompressible fluid. Thus, in order to preserve the contribution of the buoyancy force, one cannot set neither $\mathrm{Fr}=1$ nor $\mathrm{Fr}\sim 𝒪\left(1\right)$. Indeed, a possibility is that the Oberbeck–Boussinesq approximate governing Equations (1)–(3) are a limiting case of the local balance equations for a fully compressible (variable density) flow when

(39)$\begin{array}{c}\delta \to 0,\hspace{1em}\mathrm{Fr}\to 0\hspace{1em}\mathrm{with}\hspace{1em}\frac{\delta }{\mathrm{Fr}}\sim 𝒪\left(1\right).\end{array}$

The double limit defined by Equation (39) is stated, although in slightly different terms, by Hills and Roberts [6] and reported also by other authors [3,22]. It is also implicitly employed in the analysis carried out by Rajagopal, Ruzicka, and Srinivasa [2].

In the case of natural convection, where buoyancy alone causes the flow, a typical choice of $\mathrm{∆}{P}_r$ is

(40)$\begin{array}{c}\mathrm{∆}{P}_r=\frac{\mu \alpha }{L^2}.\end{array}$

Thus, the ratio $\delta /\mathrm{Fr}$ coincides with the Rayleigh number,

(41)$\begin{array}{c}\mathrm{Ra}=\frac{\delta }{\mathrm{Fr}}=\frac{{\rho }_{0}g\beta \mathrm{∆}{T}_r{L}^3}{\mu \alpha }.\end{array}$

The approach based on the limit defined by Equation (39) leads to a singular behaviour of ${\mathbf{F}}^{*}$, as one can immediately infer from Equation (35), since

$-\frac{1}{\mathrm{Fr}}({\mathbf{\nabla }}^{*}{P}^{*}+{\hat{\mathbf{e}}}_y)$

tends to infinity in the limit. Stated in these terms, the singular nature of the Oberbeck–Boussinesq limit of ${\mathbf{F}}^{*}$ seems unavoidable, if one aims to keep a finite contribution of the buoyancy force. Actually, there is just one way to avoid a singular behaviour for ${\mathbf{F}}^{*}$ in the limit defined by Equation (39). One should constrain ∇P to be equal to $-{\rho }_{0}g{\hat{\mathbf{e}}}_y$ whenever the Oberbeck–Boussinesq approximation is used. More precisely, one must have

(42)$\begin{array}{c}-{\mathbf{\nabla }}^{*}{P}^{*}-{\hat{\mathbf{e}}}_y\sim 𝒪\left(\mathrm{Fr}\right),\end{array}$

or, equivalently, by employing Equation (39),

(43)$\begin{array}{c}-{\mathbf{\nabla }}^{*}{P}^{*}-{\hat{\mathbf{e}}}_y\sim 𝒪\left(\delta \right).\end{array}$

This result can be reformulated by defining a dimensionless scalar field, ${\mathcal{P}}^{*}$, such that, when the limit (39) is to be taken, one may write

(44)$\begin{array}{c}-{\mathbf{\nabla }}^{*}{P}^{*}-{\hat{\mathbf{e}}}_y\approx -\mathrm{Fr}{\mathbf{\nabla }}^{*}{\mathcal{P}}^{*}.\end{array}$

Thus, Equation (35) yields

(45)$\begin{array}{c}{\mathbf{F}}^{*}\approx -{\mathbf{\nabla }}^{*}{\mathcal{P}}^{*}+\frac{\delta }{\mathrm{Fr}}{T}^{*}{\hat{\mathbf{e}}}_y,\end{array}$

so that the limit (39) leads to a finite result. It must be pointed out that, from a physical viewpoint, the meaning of the term ${\mathbf{\nabla }}^{*}{\mathcal{P}}^{*}$ has no connection with the difference between the pressure gradient and the hydrostatic pressure gradient, as this difference tends to zero in the limit (39).

Interestingly enough, the derivation of the Oberbeck–Boussinesq approximation presented by Rajagopal, Ruzicka, and Srinivasa [2], recently reformulated within a simplified framework by Barletta [15], leads to the same conclusion drawn above, namely that the pressure gradient must always coincide with the hydrostatic pressure gradient. We mention that the treatment presented by Rajagopal, Ruzicka, and Srinivasa [2] is based on a dimensionless scaling of the governing balance equations where the reference length is $𝒪$ and the reference velocity is $𝒪$, where $ϵ$ is the perturbation parameter. The Oberbeck–Boussinesq approximation is defined as the asymptotic case where $ϵ\to 0$. Thus, the dimensionless scaling is singular in this limit. Incidentally, $ϵ$ is proportional to ${\delta }^{1/3}$, with $\delta$ given by Equation (37). The study reported by Rajagopal, Ruzicka, and Srinivasa [2] discusses in detail the serious drawbacks of the previous theoretical studies that define the Oberbeck–Boussinesq approximation as a limiting case obtained by letting one or more perturbation parameters to zero. Examples are the papers by Spiegel and Veronis [4], by Gray and Giorgini [5] and by Hills and Roberts [6]. As a consequence, one can say that an asymptotic theory, based on a suitable perturbation scheme, which is aimed at a completely rigorous and physically convincing deduction of the Oberbeck–Boussinesq set of governing Equations (1)–(3) is based on the ideas conveyed in the paper by Rajagopal, Ruzicka, and Srinivasa [2] under the restatement provided by Barletta [15].

7. Conclusions

The Oberbeck–Boussinesq approximation for the local balance equations of mass, momentum and energy of buoyancy–induced fluid flows has been briefly outlined. It has been stressed that the linear temperature–dependent density expression is to be used, within the local momentum balance equation, only to transform the combined pressure gradient force and gravity force into a combination of terms involving the dynamic pressure gradient and the buoyancy force. Out of this very specific use, the variable density expression cannot be employed consistently with the approximated Oberbeck–Boussinesq scheme. Three examples have been discussed, relative to either natural or mixed convection flows, where it has been shown that an improper use of the variable density unavoidably leads to violations of the principle of mass conservation. A final Section has been devoted to the idea of the Oberbeck–Boussinesq approximation as a limiting case of the fully compressible local balance equations where the dimensionless parameter $\delta =\beta \mathrm{∆}{T}_r$ tends to zero, with $\beta$ the thermal expansion coefficient and $\mathrm{∆}{T}_r$ the reference temperature difference. It has been pointed out that such an idea leads either to a singular behaviour for the gradient of the dynamic pressure, or to the constraint that such a gradient be identically zero. We conclude that the Oberbeck–Boussinesq approximation is to be founded on the scenario where the dynamic pressure gradient is considered as vanishingly small, albeit of the same order of magnitude of the buoyancy force as pointed out by Barletta [15].

Author Contributions

A.B., M.C. and D.A.S.R. contributed equally to this work. All authors have read and agreed to the published version of the manuscript.

Data Availability Statement

Not applicable.

Conflicts of Interest

The authors declare no conflict of interest.

Footnotes

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Figures

Figure 1. Sketch of the rectangular cavity with side heating discussed in Section 3.

Figure 2. Sketch of the mixed convection duct flow discussed in Section 4.

Figure 3. Sketch of the porous slab separating two fluid reservoirs discussed in Section 5.

Appendix A. Batchelor’s Profile

If one considers the rectangular cavity displayed in Figure 1 and assumes $H\gg L$, namely an infinitely tall cavity, then the velocity distribution at the midplane $y=0$ can be expressed as the cubic Batchelor’s profile [23], (A1)$$\begin{array}{c}v=-\frac{{\rho }_0\beta g(T_1-T_2)}{12\mu L}x(L^2-x^2),\end{array}$$ with the linear temperature profile, (A2)$$\begin{array}{c}T-T_0=-\frac{T_1-T_2}{2L}x.\end{array}$$

Equation (A2) defines the pure conduction regime happening in the cavity. On the basis of Equation (14), Equations (A1) and (A2) yield the mass flow rate per unit depth (in the z-direction), (A3)$$\begin{array}{c}\dot{m}=-{\rho }_0\beta \underset{-L}{\overset{L}{\int }}(T-T_0)v\mathrm{d}x=-\frac{{\rho }_0^2{\beta }^{2}g{(T_1-T_2)}^2{L}^3}{90\mu }.\end{array}$$

References

1. Zeytounian, R.K. Asymptotic Modeling of Atmospheric Flows; Springer: Berlin/Heidelberg, Germany, 1990; [DOI: https://dx.doi.org/10.1007/978-3-642-73800-5]

2. Rajagopal, K.R.; Ruzicka, M.; Srinivasa, A.R. On the Oberbeck–Boussinesq approximation. Math. Model. Methods Appl. Sci.; 1996; 6, pp. 1157-1167. [DOI: https://dx.doi.org/10.1142/S0218202596000481]

3. Zeytounian, R.K. Joseph Boussinesq and his approximation: A contemporary view. C. R. Méc.; 2003; 331, pp. 575-586. [DOI: https://dx.doi.org/10.1016/S1631-0721(03)00120-7]

4. Spiegel, E.A.; Veronis, G. On the Boussinesq approximation for a compressible fluid. Astrophys. J.; 1960; 131, pp. 442-447. [DOI: https://dx.doi.org/10.1086/146849]

5. Gray, D.D.; Giorgini, A. The validity of the Boussinesq approximation for liquids and gases. Int. J. Heat Mass Transf.; 1976; 19, pp. 545-551. [DOI: https://dx.doi.org/10.1016/0017-9310(76)90168-X]

6. Hills, R.N.; Roberts, P.H. On the Motion of a Fluid That Is Incompressible in a Generalized Sense and Its Relationship to the Boussinesq Approximation. Magnetohydrodynamics and the Earth’s Core: Selected Works by Paul Roberts; Soward, A.M. CRC Press/Taylor & Francis Group: London, UK, 2003; pp. 264-271. [DOI: https://dx.doi.org/10.1201/b12581]

7. Vadasz, P. The effect of thermal expansion on porous media convection part 1: Thermal expansion solution. Transp. Por. Med.; 2001; 44, pp. 421-443. [DOI: https://dx.doi.org/10.1023/A:1010708823008]

8. Rajagopal, K.R.; Saccomandi, G.; Vergori, L. On the Oberbeck–Boussinesq approximation for fluids with pressure dependent viscosities. Nonlinear Anal. Real World Appl.; 2009; 10, pp. 1139-1150. [DOI: https://dx.doi.org/10.1016/j.nonrwa.2007.12.003]

9. Rajagopal, K.R.; Saccomandi, G.; Vergori, L. On the approximation of isochoric motions of fluids under different flow conditions. Proc. R. Soc. A Math. Phys. Eng. Sci.; 2015; 471, 20150159. [DOI: https://dx.doi.org/10.1098/rspa.2015.0159]

10. Mayeli, P.; Sheard, G.J. Buoyancy-driven flows beyond the Boussinesq approximation: A brief review. Int. Commun. Heat Mass Transf.; 2021; 125, 105316. [DOI: https://dx.doi.org/10.1016/j.icheatmasstransfer.2021.105316]

11. Arnone, G.; Capone, F.; De Luca, R.; Massa, G. Compressibility effect on Darcy porous convection. arXiv; 2022; [DOI: https://dx.doi.org/10.48550/arXiv.2210.08075]

12. Bella, P.; Feireisl, E.; Oschmann, F. The incompressible limit for the Rayleigh–Bénard convection problem. arXiv; 2022; [DOI: https://dx.doi.org/10.48550/arXiv.2206.14041]

13. Gebhart, B.; Jaluria, Y.; Mahajan, R.L.; Sammakia, B. Buoyancy–Induced Flows and Transport; Hemisphere Publishing: New York, NY, USA, 1988.

14. Straughan, B. The Energy Method, Stability, and Nonlinear Convection; Springer Science+Business Media: New York, NY, USA, 1992; Chapter 3 [DOI: https://dx.doi.org/10.1007/978-1-4757-2194-2_3]

15. Barletta, A. The Boussinesq approximation for buoyant flows. Mech. Res. Commun.; 2022; 124, 103939. [DOI: https://dx.doi.org/10.1016/j.mechrescom.2022.103939]

16. de Vahl Davis, G. Laminar natural convection in an enclosed rectangular cavity. Int. J. Heat Mass Transf.; 1968; 11, pp. 1675-1693. [DOI: https://dx.doi.org/10.1016/0017-9310(68)90047-1]

17. Barletta, A. A proof that convection in a porous vertical slab may be unstable. J. Fluid Mech.; 2015; 770, pp. 273-288. [DOI: https://dx.doi.org/10.1017/jfm.2015.154]

18. Barletta, A. Instability of stationary two–dimensional mixed convection across a vertical porous layer. Phys. Fluids; 2016; 28, 014101. [DOI: https://dx.doi.org/10.1063/1.4939287]

19. Barletta, A.; Celli, M. The Horton-Rogers-Lapwood problem for an inclined porous layer with permeable boundaries. Proc. R. Soc. A Math. Phys. Eng. Sci.; 2018; 474, 20180021. [DOI: https://dx.doi.org/10.1098/rspa.2018.0021] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/30333696]

20. Barletta, A.; Rees, D.A.S. On the onset of convection in a highly permeable vertical porous layer with open boundaries. Phys. Fluids; 2019; 31, 074106. [DOI: https://dx.doi.org/10.1063/1.5110484]

21. Vynnycky, M.; Mitchell, S.L. On steady natural convection in a vertical porous slab with permeable boundaries. J. Fluid Mech.; 2022; 935, A21. [DOI: https://dx.doi.org/10.1017/jfm.2022.16]

22. Zeytounian, R.K. The Bénard problem for deep convection: Rigorous derivation of approximate equations. Int. J. Eng. Sci.; 1989; 27, pp. 1361-1366. [DOI: https://dx.doi.org/10.1016/0020-7225(89)90060-8]

23. Batchelor, G.K. Heat transfer by free convection across a closed cavity between vertical boundaries at different temperatures. Quart. Appl. Math.; 1954; 12, pp. 209-233. [DOI: https://dx.doi.org/10.1090/qam/64563]