3.1. Governing equations

We consider the 2-D, steady flow in a $(x,z)$ cross-section of an infinitely wide duct, where $x$ is the along-channel coordinate and $z$ is the coordinate going from the bottom to the top of the duct. The fluid velocity is $\boldsymbol {v}=(u,0,w)$. The flow is described by the continuity and the steady Navier–Stokes equations with the Boussinesq approximation for an incompressible fluid,

(3.1)$$\begin{gather} \frac{\partial u}{\partial x}+\frac{\partial w}{\partial z} = 0, \end{gather}$$

(3.2)$$\begin{gather}u\frac{\partial u}{\partial x}+w\frac{\partial u}{\partial z} ={-}\frac{1}{\bar{\rho}}\frac{\partial p}{\partial x}+\nu \left(\frac{\partial^{2} u}{\partial x^{2}}+\frac{\partial^{2} u}{\partial z^{2}}\right)+g \frac{\rho'}{\bar{\rho}}\sin \theta, \end{gather}$$

(3.3)$$\begin{gather}u\frac{\partial w}{\partial x}+w\frac{\partial w}{\partial z} ={-}\frac{1}{\bar{\rho}}\frac{\partial p}{\partial z}+\nu \left(\frac{\partial^{2} w}{\partial x^{2}}+\frac{\partial^{2} w}{\partial z^{2}} \right) - g \frac{\rho'}{\bar{\rho}}\cos \theta, \end{gather}$$

where $p$ is the pressure and $\rho '=\rho -\bar {\rho }$ is the variable part of the total density, with $\bar {\rho } \gg \rho '$. In the reference frame of the inclined duct, the gravity vector is $\boldsymbol {g}=(g \sin \theta, 0, -g \cos \theta )$. A linear equation of state relates the density $\rho$ to, for example, salt concentration or fluid temperature. In this way, the density $\rho$ is also governed by a steady transport equation

(3.4)\begin{equation} u\frac{\partial \rho}{\partial x}+ w\frac{\partial \rho}{\partial z}= \kappa \left( \frac{\partial^{2} \rho}{\partial x^{2}} + \frac{\partial^{2} \rho}{\partial z^{2}} \right). \end{equation}

We define the non-dimensional variables denoted by a tilde such that

(3.5a–f)\begin{equation} u=U \tilde{u},\quad w=\frac{U}{A} \tilde{w},\quad x=\frac{L}{2}\tilde{x},\quad z=\frac{H}{2}\tilde{z},\quad \rho'={\rm \Delta}\rho \tilde{\rho}',\quad p=\frac{2\bar{\rho} U \nu L}{H^2} \tilde{p}. \end{equation}

Using the non-dimensional variables defined in (3.5a–f), (3.1)–(3.4) can be written as

(3.6)$$\begin{gather} \frac{\partial \tilde{u}}{\partial \tilde{x}}+\frac{\partial \tilde{w}}{\partial \tilde{z}}=0, \end{gather}$$

(3.7)$$\begin{gather}A^{{-}1} {{Re}}_g\, {{Fr}} \left(\tilde{u}\frac{\partial \tilde{u}}{\partial \tilde{x}}+\tilde{w} \frac{\partial \tilde{u}}{\partial \tilde{z}}\right) ={-}\frac{\partial \tilde{p}}{\partial \tilde{x}}+A^{{-}2} \frac{\partial^{2} \tilde{u}}{\partial \tilde{x}^2}+ \frac{\partial^{2} \tilde{u}}{\partial \tilde{z}^2}+ \frac{{{Re}}_g}{2 {{Fr}}} \tilde{\rho}'\sin \theta, \end{gather}$$

(3.8)$$\begin{gather}A^{{-}3} {{Re}}_g {{Fr}} \left(\tilde{u}\frac{\partial \tilde{w}}{\partial \tilde{x}}+\tilde{w} \frac{\partial \tilde{w}}{\partial \tilde{z}}\right) ={-}\frac{\partial \tilde{p}}{\partial \tilde{z}}+A^{{-}4} \frac{\partial^{2} \tilde{w}}{\partial \tilde{x}^{2}}+A^{{-}2} \frac{\partial^{2} \tilde{w}}{\partial \tilde{z}^{2}} - A^{{-}1}\frac{{{Re}}_g}{2{{Fr}}} \tilde{\rho}' \cos \theta, \end{gather}$$

(3.9)$$\begin{gather}A^{{-}1}{{Re}}_g\,{{Sc}}\,{{Fr}} \left(\tilde{u}\frac{\partial \tilde{\rho}'}{\partial \tilde{x}}+ \tilde{w}\frac{\partial \tilde{\rho}'}{\partial \tilde{z}}\right)= A^{{-}2}\frac{\partial^{2} \tilde{\rho}'}{\partial \tilde{x}^{2}} + \frac{\partial^{2} \tilde{\rho}'}{\partial \tilde{z}^{2}}, \end{gather}$$

where the Froude number is defined as

(3.10)\begin{equation} {{Fr}}\equiv \frac{U}{U_g} = \frac{1}{HU_g}\int_{{-}H/2}^{H/2}|u|\,{\rm d} z. \end{equation}

We do not set, a priori, $U_g$ as the typical velocity scale (as done by, for example, Cormack et al. Reference Cormack, Leal and Imberger1974) because we do not know if $U$ scales with $U_g$, and it is imperative to use a representative velocity scale to fully profit from the asymptotic analysis. In fact, determining the relationship between $U$ and $U_g$ (i.e. determining $Fr$) as a function of the control parameters of the problem is a crucial step in deriving the analytical solution in the HGV-A approximation. Hence, the Froude number represents the non-dimensional volume flow rate through the duct, and it is related to the function $f_{{\rm \Delta} U}$ discussed by Lefauve & Linden (Reference Lefauve and Linden2020a) for which they used the peak-to-peak velocity as the typical velocity scale.

3.2. Asymptotic analysis for long ducts

We now make use of the fact that $A^{-1}\ll 1$ to describe the main characteristics and balances in the HGV-A approximation using an asymptotic analysis. We introduce the symbol $O$ to denote the mathematical order of a function and $O_s$ to denote the sharp order (see Appendix A for formal definitions). Note that the symbol $O$ does not denote the physical order of magnitude because no account is kept of constants of proportionality (Van Dyke Reference Van Dyke1975). Our approach differs in three ways from that of Cormack et al. (Reference Cormack, Leal and Imberger1974). First, as already mentioned, we do not consider a priori $U_g$ as a representative velocity scale. Second, we simplify the problem by considering only the lowest order terms in the expansion. Third, we complicate the problem by not assuming that the interface is parallel to the bottom and top of the duct.

The starting point of the asymptotic analysis is to write the velocity, the pressure and the variable part of the density as asymptotic expansions on $A^{-1}$, such that

(3.11a–d)\begin{equation} \tilde{u}=\sum_{n=0}^\infty A^{{-}n}\tilde{u}_n, \quad \tilde{w}=\sum_{n=0}^\infty A^{{-}n}\tilde{w}_n, \quad \tilde{p}=\sum_{n=0}^\infty A^{{-}n}\tilde{p}_n, \quad \tilde{\rho}'=\sum_{n=0}^\infty A^{{-}n}\tilde{\rho}'_n. \end{equation}

As a second step, it is necessary to express the response parameters, in this case, $Fr$ as a function of $A^{-1}$. We write

(3.12)\begin{equation} {{Fr}} = C A^{{-}l}, \end{equation}

where $C$ is a finite constant such that ${{Fr}} =O_s(A^{-l})$ with $l \in \mathbb {N}$. We now substitute (3.11a–d) and (3.12) into (3.7)–(3.9), and extract the equations for the lowest (zeroth-order) terms.

We start the analysis with the $x$ and $z$ components of the Navier–Stokes equations, (3.7) and (3.8), respectively. Considering $A^{-1}\ll 1$ immediately yields, as expected for long ducts, that horizontal viscous diffusion (terms with $\partial ^2/\partial x^2$) can be neglected from these equations when compared with vertical viscous diffusion (terms with $\partial ^2/\partial z^2$) because of the additional $A^{-2}$ factor. The lowest order term of vertical viscous diffusion in (3.8) is $O(A^{-2})$ and can also be immediately neglected with respect to the lowest order terms of the pressure gradient and gravity forces that are $O(A^{0})$.

Different balances can be obtained by considering different values of $l$ in (3.12). However, here we are interested only in the laminar flow for which we assume that the nonlinear terms in (3.7) and (3.8) can be neglected, and we will verify this a posteriori by comparing with laboratory experiments. Neglecting the nonlinear terms formally means that we consider here only the case where ${{Re}}_g {{Fr}} A^{-1}=O(A^{-1})$. Equation (3.8) yields then, at lowest order, the hydrostatic balance

(3.13)\begin{equation} \frac{\partial \tilde{p}_0}{\partial \tilde{z}} ={-}\frac{K}{2} \tilde{\rho}'_0\cos \theta, \end{equation}

with

(3.14)\begin{equation} K \equiv \dfrac{{{Re}}_g}{A\,{{Fr}}}. \end{equation}

However, the hydrostatic balance (3.13) can only hold if $K\cos \theta$ is a finite constant (formally, if $K\cos \theta =O_s(A^0)$), so that the terms on both sides are of the same order. Furthermore, $\theta \ll 1$ so $\cos \theta \approx 1$, and thus, $K$ must be a finite constant ($0< K<\infty$) in the limit $A^{-1}\to 0$ (formally, $K=O_s(A^0$)).

We now turn our attention to (3.7), where the remaining terms yield

(3.15)\begin{equation} \frac{\partial^{2} \tilde{u}_0}{\partial \tilde{z}^{2}}= \frac{\partial \tilde{p}_0}{\partial \tilde{x}}-\frac{K}{2} \tilde{\rho}'_0 A \sin \theta. \end{equation}

Here, we can clearly see that the hydrostatic pressure and gravity are balanced by viscous momentum diffusion. This is the balance that gives rise to the HGV acronym. Here, the condition that $\theta \ll 1$ is also important so that $A \sin \theta$ remains finite when $A\gg 1$. Note that taking $\sin \theta =0$ yields the governing equation for a horizontal duct or what Lefauve & Linden (Reference Lefauve and Linden2020a) call the hydrostatic/viscous balance. Finally, the $z$ component of the velocity, $\tilde {w}$, is given by the continuity equation

(3.16)\begin{equation} \frac{\partial \tilde{u}_0}{\partial \tilde{x}}+ \frac{\partial \tilde{w}_0}{\partial \tilde{z}}=0. \end{equation}

The simplified equations (3.13), (3.15) and (3.16) were already introduced by Macagno & Rouse (Reference Macagno and Rouse1961) as the governing equations for SID flows. However, the perturbation analysis yields the additional important result that $K=O_s(A^0)$. This means that ${Fr=O_s(A^{-1})}$, i.e. that ${{Fr}}=U/U_g = {{Re}}_g /( A\,K)$ (in agreement with the numerical results of Kaptein et al. Reference Kaptein, van de Wal, Kamp, Armenio, Clercx and Duran-Matute2020), with $K^{-1}$ a proportionality constant in the limit of $A^{-1}\to 0$. We can see from (3.15) that the balance and, hence, the value of $K$ still depends on the value of $A \sin \theta$ that remains a finite parameter in the simplified set of equations. In short, the perturbation analysis so far tells us that, to find the dimensionless volumetric flow rate $Fr= {{Re}}_g /(K\, A)$, we must find the value of $K$ as a function of $A \sin \theta$ since the values of $Re_g$ and $A$ are known.

We now continue our analysis with the density transport equation (3.9), which yields at lowest order

(3.17)\begin{equation} \frac{\partial^{2} \tilde{\rho}'_0}{\partial \tilde{z}^{2}} = 0, \end{equation}

meaning that the vertical diffusion of density is equal to zero.

In short, the HGV-A approximation is characterized by a hydrostatic balance in the vertical, while in the horizontal, the flow is driven by a pressure gradient force, and by gravity in the along-duct direction if $\theta \neq 0$. Finally, we have the hallmark of the HGV-A approximation: diffusion dominates over inertia in the along-channel momentum equation but is equal to zero in the density transport equation. This combination requires that $Sc$ is large (as observed by Kaptein et al. Reference Kaptein, van de Wal, Kamp, Armenio, Clercx and Duran-Matute2020), which is the case, for example, for water flows where the density difference is caused by a difference in salt concentration (${{Sc}}\approx 700$) but not by a difference in temperature ($Sc\approx 7$).

3.3. Analytical solution in the HGV-A approximation

We consider now the governing equations at lowest order (3.13) and (3.15)–(3.17), referred to as the HGV-A approximation. For simplicity in the notation, we drop the tildes. The goal now is to find an analytical solution for this set of equations. We split the procedure in two steps. First, we determine the expressions for the along-channel velocity $u_0$ and the density $\rho '_0$. Second, we determine the value of $Fr$ that, as mentioned earlier, is equivalent to finding the value of $K$ as a function of $A \sin \theta$.

3.3.1. Determining the along-channel velocity and the density

To find an expression for the along-channel velocity $u_0$ and the density $\rho '_0$ in the HGV-A approximation, we begin with the transport equation (3.17). For the solution to this equation, we propose

(3.18)\begin{equation} \rho_0'(x,z) = \frac{1}{2} - \mathcal{H}(z-\eta(x)), \end{equation}

where $\mathcal {H}(z)$ is the Heaviside function defined as $\mathcal {H}(z>0)=1$, $\mathcal {H}(z<0)=0$ and $\mathcal {H}(z=0)=1/2$. The physical meaning of this solution is that the density gets organized in two layers with a sharp interface located at $z=\eta (x)$ with $\eta (0)=0$ due to the definition of the coordinate system.

The proposed solution satisfies $\partial ^2 \rho '_0/\partial z^2=0$ in each of the layers. However, formally, the solution is not defined at the interface, and it is convenient to approach the problem by treating each layer separately, with the interface as a boundary between them. We use the index $\zeta =\pm 1$ with $\zeta =+1$ referring to the top layer and $\zeta =-1$ to the bottom layer. Using this notation, $\rho '_{0,\zeta }=-\zeta /2$. Macagno & Rouse (Reference Macagno and Rouse1961) also solved the same governing equations for a two-layer system, but they considered the simplified problem where the interface is parallel to the top and bottom of the duct, i.e. $\eta (x)=0$.

We now turn our attention to (3.13) to derive the pressure distribution and the horizontal pressure gradient that drives the flow. Integrating this equation with respect to $z$ yields an expression for the pressure,

(3.19)\begin{equation} p_{0,\zeta}(x,z)=K\frac{\zeta}{4}z\cos \theta + \gamma_\zeta(x), \end{equation}

with $\gamma _\zeta (x)$ as integration constants that are functions of $x$ to be determined using the boundary conditions. At $z=1$, the pressure is an unknown function of $x$, but it is convenient to define, without loss of generality,

(3.20)\begin{equation} p_{0,+1} (x,z=1)={-}K[f(x)- \eta(x) \cos\theta]/4, \end{equation}

with $f(x)$ an unknown function giving the $x$ dependence of the barotropic pressure. Using (3.19) and (3.20), one obtains an expression for $\gamma _{+1}(x)$ as a function of $\eta (x)$ and $f(x)$, and applying continuity of pressure at $z=\eta (x)$ yields an expression for $\gamma _{-1}(x)$. The pressure as a function of $x$ and $z$ is then given by

(3.21)\begin{equation} p_{0,\zeta}(x,z)=\zeta\frac{K}{4}[\cos\theta(z-\eta(x)-\zeta)+\zeta f(x)] \end{equation}

so that the along-channel pressure gradient is

(3.22)\begin{equation} \frac{\partial p_{0,\zeta}}{\partial x}=\zeta\frac{K}{4}\left(\cos\theta \frac{{\rm d} \eta(x)}{{\rm d}\kern0.06em x}-\zeta \frac{{\rm d} f(x)}{{\rm d}\kern0.06em x}\right), \end{equation}

where we can see that it is composed of a baroclinic part due to the sloping interface and a barotropic part given by the gradient of $f(x)$.

Substituting (3.22) into (3.15) yields

(3.23)\begin{equation} \frac{\partial^2 u_{0,\zeta}}{\partial z^2}={-}\zeta\frac{K}{4} F_\zeta(x), \end{equation}

with

(3.24)\begin{equation} F_\zeta(x)=\cos\theta \frac{{\rm d} \eta(x)}{{\rm d}\kern0.06em x}-\zeta\frac{{\rm d} f(x)}{{\rm d}\kern0.06em x}-A\sin \theta. \end{equation}

Integrating twice with respect to $z$ and applying the boundary conditions $u_0(z=\pm 1) =u_0(z=\eta (x))=0$ gives

(3.25)\begin{equation} u_{0,\zeta}(x,z) ={-}\zeta\frac{K}{8}F_\zeta(x)(z-\zeta)(z-\eta(x)). \end{equation}

Since the barotropic pressure gradient is responsible for an equal volume flow rate in both directions, we use the condition of zero mean flow through the duct,

(3.26)\begin{equation} \int_{{-}1}^{1} u_0(x,z)\,{\rm d} z=0, \end{equation}

to determine ${\rm d} f(x)/{{\rm d}\kern0.06em x}$. We do this by integrating $u_{0,-1}(x,z)$ from $z=-1$ to $z=\eta (x)$ and $u_{0,+1}(x,z)$ from $z=\eta (x)$ to $z=1$, yielding

(3.27)\begin{equation} \frac{{\rm d} f(x)}{{\rm d}\kern0.06em x}={-}\frac{\eta(x)[3+\eta^2(x)]}{1+3\eta^2(x)} \left(\cos\theta \frac{{\rm d} \eta(x)}{{\rm d}\kern0.06em x}-A\sin \theta \right), \end{equation}

and $F_\zeta (x)$ can be rewritten as

(3.28)\begin{equation} F_\zeta(x)= \frac{[1+\zeta\eta(x)]^3}{1+3\eta^2(x)} \left(\cos \theta \frac{{\rm d} \eta(x)}{{\rm d}\kern0.06em x}-A\sin \theta\right). \end{equation}

Replacing this expression for $F_\zeta$ into (3.23) clearly shows the physical origin of the two drivers of the flow: (i) a baroclinic pressure gradient due to a sloping interface between the two layers with homogeneous density, and (ii) gravity due to the tilt of the duct. For a given value of $A$, the relative importance of these two forcing terms varies with the angle $\theta$ and the slope of the interface. For example, for horizontal ducts, the pressure gradient will be the only forcing, while it would be expected that gravity takes over with increasing $\theta$ values (particularly, if ${\rm d} \eta (x)/{{\rm d}\kern0.06em x}$ decreases simultaneously). Hence, considering these two forcing terms without neglecting a priori any of them is crucial to link our knowledge of horizontal and inclined ducts. The difficulty here is that although the value of $\theta$ is known since it is a control parameter, the value of ${\rm d} \eta (x)/{{\rm d}\kern0.06em x}$ is not.

Now, the only missing part of the solution is to determine the shape of the interface $\eta (x)$. For this, we consider the flow rate through each of the layers yielding

(3.29)\begin{equation} -\int_{\eta(x)}^{1} u_{0,+1}(x,z)\,{\rm d} z=1 =\int_{{-}1}^{\eta(x)}u_{0,-1}(x,z)\,{\rm d} z, \end{equation}

where the value of one is due to the way the velocity was made dimensionless. In this way, we obtain an autonomous differential equation for $\eta (x)$,

(3.30)\begin{equation} \cos\theta \frac{{\rm d} \eta(x)}{{\rm d}\kern0.06em x}={-}\frac{48}{K} \frac{1+3\eta^2(x)}{[1-\eta^2(x)]^3}+A\sin \theta. \end{equation}

Equation (3.30) can be further used to rewrite $F_\zeta (x)$ in (3.28) yielding

(3.31)\begin{equation} F_\zeta(x)={-}\frac{48}{K}\frac{1}{[1-\zeta\eta(x)]^3}, \end{equation}

such that (3.23) is now written as

(3.32)\begin{equation} \frac{\partial ^2 u_{0,\zeta}}{\partial z^2}={-}\frac{12\zeta}{[1-\zeta\eta(x)]^3}. \end{equation}

Finally, the along-duct velocity given by (3.25) can be written as

(3.33)\begin{equation} u_{0,\zeta}(x,z) =\frac{6 \zeta}{[1-\zeta\eta(x)]^3}(z-\zeta)(z-\eta(x)), \end{equation}

where we can notice that the velocity profile at $x=0$,

(3.34)\begin{equation} u_{0,\zeta}(0,z) =6 \zeta (z-\zeta)z, \end{equation}

always has the same shape consisting of two parabolas. To reach the complete solution, it is still necessary to determine the position of the interface, which is equivalent to finding the value of $K$ as a function of $A \sin \theta$, which can be done by solving (3.30) for $\eta (x)$.

3.3.2. Determining $K$

Since the analytical solution for horizontal ducts is tractable, we focus first on this case for which (3.30) simplifies to

(3.35)\begin{equation} \frac{{\rm d} \eta(x)}{{\rm d}\kern0.06em x}={-}\frac{48}{K}\frac{1+3\eta(x)^2}{[1-\eta(x)^2]^3}, \end{equation}

yielding

(3.36)\begin{equation} x={-}\frac{K}{19\,440}[320\sqrt{3} \arctan (\sqrt{3} \eta(x))-555 \eta(x)+150 \eta(x)^3 - 27 \eta(x)^5]. \end{equation}

The value of $K$ is obtained by imposing boundary conditions. Finding the slope of the interface from an autonomous equation – similar to (3.30) – was previously done for horizontal ducts by Gu & Lawrence (Reference Gu and Lawrence2005) and for inclined ducts by Lefauve & Linden (Reference Lefauve and Linden2020a) who derived the equation using internal hydraulics. Gu & Lawrence (Reference Gu and Lawrence2005) found in fact a similar expression to (3.36), but the constant to be determined from the boundary conditions was the composite Froude number given a certain magnitude of the frictional effects. When using internal hydraulics, the boundary condition is proposed using maximum exchange flow theory (Armi & Farmer Reference Armi and Farmer1986). This theory states that the maximum flow rate possible is such that the flow becomes critical at the ends of the channel (i.e. that the composite Froude number becomes unity at $x=\pm 1$; see e.g. Dalziel Reference Dalziel1991; Zaremba, Lawrence & Pieters Reference Zaremba, Lawrence and Pieters2003; Gu & Lawrence Reference Gu and Lawrence2005 for details). In such a case, the flow is said to be hydraulically controlled. However, enforcing the composite Froude number to be equal to unity at the edges is inconsistent with the results of the perturbation analysis leading to the HGV-A approximation. On the one side, for internal hydraulics to be valid, inertia in the $x$ component of the momentum equation should not be negligible. On the other hand, setting the value of the composite Froude number at $x=\pm 1$ would mean that the position of the interface at these locations depends on the Froude number. Imposing this boundary condition to find the value of $K$ would, in turn, make $K$ depend on the Froude number. However, this is inconsistent with the results from the asymptotic analysis that $K=O_s(A^0)$ while ${{Fr}}=O_s(A^{-1})$, which is a condition needed for the hydrostatic balance (3.13) to hold.

We propose then a different boundary condition inspired by the results of Kaptein et al. (Reference Kaptein, van de Wal, Kamp, Armenio, Clercx and Duran-Matute2020) for the simulation used to exemplify the HGV-A approximation (${{Re}}_g=500$, ${A=60}$, ${{Sc}}=300$, $\theta =0^\circ$) and shown in figure 2. In this figure we can see that the interface curves sharply when reaching the end of the duct, with the light fluid turning upwards at $x=-1$ and the dense fluid turning downwards at $x=1$. The currents at the ends of the duct must have a dimensionless thickness $\epsilon$ in the $x$ direction with $\epsilon < A^{-1}\ll 1$, meaning that they are thin with respect to the length of the duct. We propose then to impose $\eta [x=\pm ( 1+\epsilon )]=\mp 1$. For long ducts, we can assume a small error of order $A^{-1}$, and enforce, instead $\eta (x=\pm 1)=\mp 1$, meaning also that ${\rm d}\eta /{{\rm d}\kern0.06em x}|_{x=\pm 1}=\mp \infty$. Note that, even though the flow inside the duct is not controlled by the composite Froude number being equal to one at the edges, the boundary conditions at the edges remain crucial in determining the slope of the interface, and hence, the volume flow rate. Taking $\eta (\pm 1)=\mp 1$ yields $K\approx 131$ and ${\rm d}\eta /{{\rm d}\kern0.06em x} |_{x=0}=-48/K \approx -0.366$ for a horizontal duct. Note that the value of ${\rm d}\eta /{{\rm d}\kern0.06em x} |_{x=0}$ is close to the slope of $-1/3$ obtained empirically and numerically by Kaptein et al. (Reference Kaptein, van de Wal, Kamp, Armenio, Clercx and Duran-Matute2020) and shows good agreement with the density field shown in figure 2.

Figure 2.

Density field from the numerical simulation by Kaptein et al. (Reference Kaptein, van de Wal, Kamp, Armenio, Clercx and Duran-Matute2020) for ${{Re}}_g=500$, $A=60$, ${{Sc}}=300$ and $\theta =0$. The white dashed lines represent the limits of the duct ($x=\pm 1$). The black dashed line represents the interface given by (3.36) with $K=131$.

To find the value of $K$ for the inclined ducts, we follow a similar approach as for the horizontal ducts, but we do the calculations numerically. The value of $K$ is obtained by solving (3.30) while enforcing $\eta (\pm 1)=\mp 1$. We further simplify the problem by taking $\cos \theta \approx 1$ since $\theta \ll 1$. It was already discussed in § 3.2 that $K$ is a constant for $A^{-1}\to 0$ but that it depends on $A \sin \theta$, and this can be seen again in (3.30). The solution to (3.30) for three different values of $A \sin \theta$ is shown in figure 3. Several properties of the shape of the interface expected from previous work (see e.g. Gu & Lawrence Reference Gu and Lawrence2005; Lefauve et al. Reference Lefauve, Partridge and Linden2019a; Kaptein et al. Reference Kaptein, van de Wal, Kamp, Armenio, Clercx and Duran-Matute2020) and observed in figure 2 are reproduced by the analytical solution. For example, the solution naturally yields a slope of the interface that is constant to a good approximation over a large portion of the duct around $x=0$, and that it bends up or down as it approaches the edges of the duct. Furthermore, the slope of the interface at $x=0$ is given by

(3.37)\begin{equation} S\equiv \left.\frac{{\rm d} \eta(x)}{{\rm d}\kern0.06em x}\right|_{x=0}\approx\cos \theta \left.\frac{{\rm d} \eta(x)}{{\rm d}\kern0.06em x}\right|_{x=0}={-}\frac{48}{K}+A \sin\theta \leq 0, \end{equation}

and it goes towards zero for increasing values of $\theta$. The values of $K$ and $S$ as a function of $A\sin \theta$ are shown in figure 4. Note that, due to the way the variables where made dimensionless, the ‘real’ slope (e.g. as observed in the experiments) is given by $S'=S/A$.

Figure 3.

Shape of the interface for three different values of $A\sin \theta$ as obtained from solving (3.30).

Figure 4.

(a) Value of the parameter $K$ defined in (3.14) and (b) value $S$, the slope of the interface at $x=0$, as a function of $A \sin \theta$. The value of $K$ and $S$ are such that the solution to the autonomous equation for $\eta (x)$ (3.30) satisfies $\eta (\pm 1)=\mp 1$ assuming $\cos \theta \approx 1$. The dotted line in (b) represents the empirical approximation used by Kaptein (Reference Kaptein2021).

Previous work has suggested a change in behaviour around $A\sin \theta =1$ where the flow transitions from lazy to forced (Lefauve et al. Reference Lefauve, Partridge and Linden2019a). In particular, the slope of the interface is considered relatively flat throughout the duct ($S\approx 0$) for $A \sin \theta >1$. Although here $S$ does not tend to zero for $A \sin \theta >1$, the values of $S'$ are quite small for typical values of $A$ used in the experiments. For example, in the case of $A=30$, the height of the interface varies about 3 mm over 1 m, which could be imperceptible by the eye. Furthermore, we can see that, for $A \sin \theta >1$, the value of $K$ varies little when compared with the variation for $A \sin \theta <1$.

The fact that we obtain the value of $S$ as a function of $A \sin \theta$ is a critical difference with respect to the parametrization used by Kaptein (Reference Kaptein2021) to determine the regime transition curves. In that work, the variation of $S$ with $A \sin \theta$ was assumed based on previous results by Lefauve et al. (Reference Lefauve, Partridge and Linden2019a) and Kaptein et al. (Reference Kaptein, van de Wal, Kamp, Armenio, Clercx and Duran-Matute2020). In particular, it was assumed that $S= (A \sin \theta -1)/3$ for $A\sin \theta <1$ and that $S=0$ for $A\sin \theta \geq 1$. Although certain features and trends of the regime transition curves might be reproduced using these assumptions, differences are also expected when the value of $S$, as shown in figure 4(b), is considered.

3.4. Implications of the HGV-A approximation for the regime transition

The derivation of the analytical solution in the HGV-A approximation yielded that $K$ is a constant for a given value of $A\sin \theta$ by solving (3.30) while enforcing $\eta (\pm 1)=\mp 1$. For the upcoming discussion and an easier comparison with the work by Lefauve et al. (Reference Lefauve, Partridge and Linden2019a), it is convenient to use the fact that, for a long duct with a small tilt angle ($\theta,A^{-1}\ll 1$), $A\sin \theta \approx \theta /\alpha$. In such a case, the solution in the HGV-A approximation yields that $K = {{Re}}_g/ ({{Fr}} A)$ must be a constant for a given value of $\theta /\alpha$. Both ${{Re}}_g$ and $A$ are control parameters of the problem, while ${{Fr}}$ is a response parameter equivalent to the non-dimensional volume flow rate. Hence, if the value of ${{Re}}_g A^{-1}$ is increased while keeping $\theta /\alpha$ fixed, the value of ${{Fr}}$ should also increase keeping the value of $K$ constant.

It is, here, convenient to define the Froude number in the HGV-A approximation as

(3.38)\begin{equation} {{Fr}}^* \equiv \frac{{{Re}}_g }{A K}\approx \frac{{{Re}}_g \alpha}{K}. \end{equation}

Note that ${{Fr}}^*$ can be seen as either a response parameter or a control parameter. It is a response parameter because it represents the non-dimensional volumetric flow rate that results from the choice of control parameters: ${{Re}}_g$, $\theta$ and $A$. In addition, it is a control parameter because ${{Re}}_g$ and $A$ are directly imposed for a given experiment, and $K$ is a geometrical parameter set by imposing the value of $\theta /\alpha$. The value of $K$ is determined from the analytical solution in the HGV-A approximation, and it is, to a very good approximation for $\theta \ll 1$, the value shown in figure 4(a). In this way, the value of ${{Fr}}^*$ is a known quantity for a given experiment. Within the HGV-A approximation, ${{Fr}} = {{Fr}}^*$, but this is not the case if the approximation does not hold. Furthermore, since $K\approx 131$ is a constant for $\theta = 0$, saying that ${{Fr}}^*$ is the control parameter for horizontal ducts is equivalent to saying that ${{Re}}_g\,A^{-1}$ is the control parameter as shown by Hogg et al. (Reference Hogg, Ivey and Winters2001).

To study the limit of validity of the HGV-A approximation, we consider a hypothetical SID set-up with given values for $A$ and $\theta$ satisfying $A\gg 1$ and $\theta \ll 1$. As just mentioned, the value of $K$ is known then. In this hypothetical set-up, we first take a sufficiently small value of $Re_g$ so that ${{Fr}}^*$ is sufficiently small to neglect inertia from the $x$ component of the momentum equation (3.7). Furthermore, $Sc\gg 1$ so that the flow gets organized (too good approximation) in a two-layer configuration and the HGV-A approximation holds. In such a case, $Fr={{Fr}}^*$. We now do a series of experiments increasing ${{Re}}_g$, which is equivalent to increasing ${{Fr}}^*$. We know from internal hydraulic theory that there is a maximum flow rate possible through the duct, i.e. a maximum possible value for $Fr$. In the frictionless case, this maximum value is $Fr=1$ that is known as the hydraulic limit (Hogg et al. Reference Hogg, Ivey and Winters2001; Lefauve et al. Reference Lefauve, Partridge and Linden2019a), while friction reduces the maximum possible value of $Fr$ (Gu & Lawrence Reference Gu and Lawrence2005). Hence, there must be a critical value ${{Fr}}^*_c$ at which the HGV-A approximation becomes invalid. For $Fr^*$ values larger than this critical value, the actual Froude number is such that $Fr< Fr^*$. Approaching this limit of validity of the HGV-A approximation, it is expected that a transition (i.e. a qualitative change in the flow) occurs.

Now, the next step is to determine the range of validity of the HGV-A approximation by determining ${{Fr}}_c^*$ since the asymptotic analysis does not give information about the physical order of magnitude at which the different terms become relevant (Van Dyke Reference Van Dyke1975). Hence, to determine the range of validity of the HGV-A approximation, it is necessary to use numerical simulations or laboratory experiments. In the following section we study this range of validity and the consequence for regime transitions in four different experimental set-ups.

4.1. Description of the data

We use the experimental data sets by Meyer & Linden (Reference Meyer and Linden2014), Lefauve et al. (Reference Lefauve, Partridge and Linden2019a) and Lefauve & Linden (Reference Lefauve and Linden2020a). These data sets have been discussed and made available online (https://doi.org/10.17863/CAM.48821, https://doi.org/10.17863/CAM.41410) by Lefauve, Partridge & Linden (Reference Lefauve, Partridge and Linden2019b) and Lefauve & Linden (Reference Lefauve and Linden2020a,Reference Lefauve and Lindenb). Four ducts with different dimensions were used. We will refer to them as mSID (m for mini), tSID (t for tall), LSID (L for large), and HSID (H for half), in accordance with Lefauve & Linden (Reference Lefauve and Linden2020a). Table 1 summarizes the characteristics of all the set-ups. For each set-up, the inclination angle of the duct $\theta$ and the gravitational Reynolds number ${{Re}}_g$ were varied independently. The fluid used was a salt (NaCl) solution ($Sc\approx 700$), and the value of ${{Re}}_g$ was varied by changing the salt concentration in each of the tanks. In total, we used 738 data points to study the regime transitions.

Table 1.

Characteristics of the experiments used in this paper. Four duct geometries (abbreviated mSID, LSID, HSID, tSID) are used (Lefauve & Linden Reference Lefauve and Linden2020a). We list the values of the dimensionless numbers describing each duct geometry ($A$ and $B$), the value of $Sc$ for salt in water and the ranges of $\theta$ and ${{Re}}_g$ explored.

Meyer & Linden (Reference Meyer and Linden2014) distinguished four different regimes: laminar (${\textsf {L}}$), Holmboe waves (${\textsf {H}}$), intermittently turbulent (${\textsf {I}}$), and turbulent (${\textsf {T}}$). Lefauve & Linden (Reference Lefauve and Linden2020a) introduced a ‘waves’ (${\textsf {W}}$) regime where waves other than Holmboe waves were observed. The different regimes were mostly identified by shadowgraph observations over a subsection of the duct, following the qualitative description of each regime by Meyer & Linden (Reference Meyer and Linden2014). A schematic of the shadowgraph set-up was presented by Lefauve (Reference Lefauve2018). The observed regime and the mass flow rate as a function of the governing parameters are provided by Lefauve & Linden (Reference Lefauve and Linden2020b). A complication for the comparison with the experimental data is that, as pointed out by Lefauve & Linden (Reference Lefauve and Linden2020a), there is a surprising difference between the results from the LSID and the mSID set-ups: the regions in the $(\theta,{{Re}}_g)$ plane where the different regimes occur do not coincide even though the values of all the dimensionless parameters are identical.

To show that the proposed solution in the HGV-A approximation is observed experimentally and to observe the limit of validity, we analyse in the following section the along-duct velocity and the density fields for three experiments in the mSID duct all with $\theta =2^\circ$ (provided by Lefauve et al. Reference Lefauve, Partridge and Linden2019b). The first experiment falls within the ${\textsf {L}}$ regime with ${{Re}}_g=398$, the second within the ${\textsf {H}}$ regime with ${{Re}}_g=1059$ and the third within the ${\textsf {I}}$ regime with ${{Re}}_g=1466$. A detailed description of the experiments and the methodology is given by Lefauve et al. (Reference Lefauve, Partridge and Linden2019a). With this choice of experiments, we have the same approach as mentioned for a hypothetical SID set-up in § 3.4, but with a real set-up and real experimental results.

4.2. Experimental confirmation of the HGV-A approximation

For the three experiments for which we analyse the velocity and density field, we first determine the value of $K$ numerically as explained in § 3.3.2. This value is given to good approximation in figure 4(a). For $\theta = 2^\circ$ and $A=30$, $A \sin \theta \approx \theta /\alpha \approx 1.0$ so $K\approx 39$. We then compute the value of the Froude number in the HGV-A approximation $Fr^*\equiv {{Re}}_g (KA)^{-1}$, giving ${{Fr}}^*\approx 0.34$ for the ${\textsf {L}}$ experiment, ${{Fr}}^*\approx 0.91$ for the ${\textsf {H}}$ experiment and ${{Fr}}^*\approx 1.25$ for the ${\textsf {I}}$ experiment. Here, we note already that the value of ${{Fr}}^*$ for the ${\textsf {I}}$ experiments is larger than ${{Fr}}=1$ corresponding to the hydraulic limit. In the following, we investigate if the key flow characteristics and working assumptions leading to the solution in the HGV-A approximation derived in § 3.3 hold in these experiments. We consider four aspects: (i) the density field, focusing on the slope of the interface; (ii) the $x$ component of the momentum equation (4.1), focusing on the fact that the HGV balance holds; (iii) the vertical profiles of the along-channel velocity component; and (iv) the prediction of the mass flux through the duct.

Figure 5 shows the experimentally obtained time-averaged density fields and the location of the interface as predicted by the HGV-A approximation for the three individual experiments considered. The slope of the interface close to $x=0$ is given to good approximation by (3.37), which yields $S\approx -0.20$ for these experiments. First, we see that the density is indeed organized into two layers with a sharp interface in between, with the ${\textsf {I}}$ experiment presenting a thicker interface. For the ${\textsf {L}}$ experiment, the semi-analytical solution for the slope of the interface agrees well with the experiment. This means that (3.30) with the boundary condition $\eta (\pm 1)=\mp 1$ predicts well the shape of the interface for this experiment and that the value of $K$ obtained from the analytical solution in the HGV-A approximation applies. The agreement is less for the ${\textsf {H}}$ experiment, with the slope of the interface being slightly less steep than predicted. Finally, for the ${\textsf {I}}$ experiment, the interface is fully parallel to the top and bottom of the duct, suggesting that the solution in the HGV-A approximation does not hold for this experiment.

Figure 5.

Time-averaged density field in the $(x,z)$ plane at $y=0$ for three experiments in the mSID set-up. (a) Experiment within the ${\textsf {L}}$ regime ($\theta =2^\circ$, ${{Re}}_g=398$). (b) Experiment within the ${\textsf {H}}$ regime ($\theta =2^\circ$, ${{Re}}_g=1059$). (c) Experiment within the ${\textsf {I}}$ regime ($\theta =2^\circ$, ${{Re}}_g=1455$). The black dashed line represents the line $z=Sx$ with $S=-0.20$ as given by (3.37).

As a next step, we investigate if the HGV balance (3.32) is indeed the leading order balance in the $x$ component of the momentum equation (3.7) close to the centre of the duct ($x\approx 0$), which we rewrite here as

(4.1)\begin{align} \underbrace{A^{{-}1} {{Re}}_g {{Fr}} \left(u\frac{\partial u}{\partial x}+w \frac{\partial u}{\partial z}\right)}_{I} = \underbrace{A^{{-}2} \frac{\partial^{2} u}{\partial x^2}}_{II}+ \underbrace{\frac{\partial^{2} u}{\partial z^2}}_{III}+ \underbrace{\frac{{{Re}}_g }{2A{{Fr}}}\left( A \sin \theta\rho'-\cos\theta \int_{0}^{z}\frac{\partial\rho'}{\partial x}{\rm d} z'\right)}_{IV}, \end{align}

with $z'$ a dummy variable. Here, we have assumed that the hydrostatic balance (3.13) holds and followed the derivation by Kaptein et al. (Reference Kaptein, van de Wal, Kamp, Armenio, Clercx and Duran-Matute2020) to write the pressure gradient term in this form. To see that the HGV balance (3.32) holds, vertical viscous diffusion (term III) must be equal to minus the forcing (term IV). In addition, we want to see that term I (inertia) is negligible with respect to term III and term IV. Figure 6 shows with black lines the vertical profiles of terms I, III and (minus) IV of (4.1) obtained from the time- and $y$-averaged velocity and density fields at $x\approx -0.2$. We do not show term II because it is practically equal to zero for all experiments. In general, we see that term I is also much smaller than the other two terms meaning that it can indeed be neglected. As can be seen, term III and (minus) term IV overlap well for all experiments, with the experiment in the ${\textsf {L}}$ regime showing better agreement, meaning that the HGV balance holds particularly well for the ${\textsf {L}}$ experiment. In the case of the ${\textsf {H}}$ and ${\textsf {I}}$ experiments, terms III and IV match in general, suggesting that the HGV balance is still the leading balance of the time-averaged flow, but some deviations are observed.

Figure 6.

The vertical profiles of the time- and $y$-averaged terms of the horizontal momentum balance (4.1) as a function of $z$ and $x\approx -0.2$ for the same three experiments as in figure 5. The black dotted, solid and dashed lines denote the terms I, III and (minus) IV, respectively, as obtained from the laboratory experiments. Term II is always approximately equal to zero and not explicitly shown. The dashed grey line represents both terms III and (minus) IV as obtained from the analytic solution in the HGV-A approximation that assumes an infinitely sharp interface, and is given by the right-hand side of (3.32).

To understand these deviations, we make a comparison with the theory. Figure 6 also shows with dashed grey lines the value for terms III and IV according to the solution in the HGV-A approximation derived in § 3.3 and given by the right-hand side of (3.32), which assumes an infinitely thin interface. The interface in the experiments is not infinitely thin due to diffusion in the vertical direction that is fully neglected in the HGV-A approximation. This effect is most visible in the ${\textsf {I}}$ experiment due to turbulent diffusion across the interface. For all three experiments, the forcing corresponds well with forcing expected in the HGV-A approximation given by the sum of gravity and the baroclinic pressure gradient. It is the diffusion in the vertical direction that shows deviations for the ${\textsf {H}}$ and ${\textsf {I}}$ experiments. We note, in particular, an increase of the magnitude of vertical momentum diffusion in the bottom boundary layer, and just above and below the shear layer. Note that the reasons why the momentum balance is not closed for the ${\textsf {H}}$ and ${\textsf {I}}$ experiments is that (4.1) would need additional terms for the time-averaged experimental results. Most importantly, it misses the momentum diffusion in the across-duct direction and a term containing the time-averaged effect of temporal perturbation (similar to Reynolds stresses).

We now compare the velocity profiles. Figure 7 shows the vertical profiles of the $x$ component of the velocity at $x\approx -0.2$ and averaged in the $y$ direction for the same three experiments. Choosing different values of $x$ does not affect the results beside the fact that the height of the point with $u=0$ is slightly shifted. For the experiment in the ${\textsf {L}}$ regime, the flow is steady so the instantaneous and time-averaged velocity profiles coincide. Even for the ${\textsf {H}}$ and ${\textsf {I}}$ experiments, the variation in the instantaneous profiles is small compared with the magnitude of the time-averaged velocity. For the three experiments, we see a different level of disagreement in the amplitude of the velocity. Computing the value of ${{Fr}}$ for the experiments using the time- and $y$-averaged $x$-velocity component yields ${{Fr}}=0.40$ (compared with ${{Fr}}^*\approx 0.34$) for the ${\textsf {L}}$ experiment, ${{Fr}}=0.81$ (compared with ${{Fr}}^*\approx 0.91$) for the experiment in the ${\textsf {H}}$ regime and ${{Fr}}=0.90$ (compared with ${{Fr}}^*\approx 1.25$) for the experiment in the ${\textsf {I}}$ regime. As mentioned earlier, maximum exchange theory says that the maximum value is ${{Fr}}=1$, so we could have already expected that the Froude number for the ${\textsf {I}}$ experiment was not going to reach ${{Fr}}^*=1.25$, but we have here confirmation of this. It suggests again that the HGV-A approximation does not hold for this experiment. On the other hand, the value of the Froude number for the ${\textsf {L}}$ and ${\textsf {H}}$ experiments are well predicted (within 13 %) by the HGV-A analytical solution. We believe that this is a very good result, particularly if we consider, on the one side, experimental errors, and on the other, some of the assumptions made during the analytical derivation. For example, we assumed an infinitely wide duct, while the experiments have a cross-sectional aspect ratio $B=1$.

Figure 7.

Vertical velocity profiles of the $x$ component of the velocity as a function of $z$ for $x\approx -0.2$ and averaged in the $y$ direction for the same three experiments as in figures 5 and 6. The grey lines represent the instantaneous velocity, the black dashed line the time-averaged velocity and the red line represents the velocity in the HGV-A approximation given by (3.33).

In the previous paragraphs, we studied the characteristics of only three experiments with equal value of $\theta /\alpha$ and, hence, equal value of $K$. In general, we saw that the ${\textsf {L}}$ experiment is very well reproduced by the HGV-A approximation, the ${\textsf {H}}$ experiment is slightly less so, and the ${\textsf {I}}$ experiment shows clear differences. In § 3.4 we discussed that we expect that, for ${{Fr}}^*$ values larger than a critical value ${{Fr}}^*_c$, the HGV-A approximation does not hold, and the actual Froude number must be smaller than ${{Fr}}^*$. The results suggest that ${{Fr}}^*_c$ is found somewhere between the ${\textsf {H}}$ and ${\textsf {I}}$ experiment. However, it is difficult to determine clear trends from only these three experiments. Although values of $Fr$ are not available for many experiments, the mass flux

(4.2)\begin{equation} Q_m\equiv\int_{{-}1}^{1}\int_{{-}1}^{1}|\rho' u|\,{\rm d} y\,{\rm d} z \end{equation}

is available for several experiments. The value of $Q_m$ is commonly used to compare (and, thus, characterize) the flow, with $Q_m=1/2$ corresponding to the hydraulic limit (Hogg et al. Reference Hogg, Ivey and Winters2001; Lefauve et al. Reference Lefauve, Partridge and Linden2019a), while $Q_m={{Fr}}/2={{Fr}}^*/2$ in the HGV-A approximation. The factor $1/2$ arises from the fact that $\rho '=\pm 1/2$ in the HGV-A approximation following the way we made the density dimensionless.

The value of $Q_m$ as a function ${{Fr}}^*/2$ for the experiments in the mSID set-up are shown in figure 8. The values of $Q_m$ can then be compared with the expected values in the HGV-A approximation ($Q_m={{Fr}}^*/2$) and to the hydraulic limit ($Q_m=1/2$). For ${{Fr}}^*/2\lesssim 1/2$, $Q_m\approx {{Fr}}^*/2$, showing that the HGV-A approximation gives a good estimate for the value of $Q_m$. Although the values for some specific experiments are overestimated, and some others are underestimated, clear trends are visible. For example, for the experiments with $\theta \geq 1.0^\circ$, $Q_m\approx {{Fr}}^*/2$ up to ${{Fr}}^*/2=1/2$. Then, a sharp transition to the hydraulic limit occurs, and $Q_m\approx 1/2$ for ${{Fr}}^*/2\geq 1/2$. For experiments with $\theta = 0.5^\circ$, $Q_m$ seems to saturate reaching a maximum value $Q_m\approx 0.4<1/2$. The results for $\theta = 0^\circ$ suggest that $Q_m$ will similarly not reach the value $1/2$ for a horizontal duct. The impossibility of reaching $Q_m=1/2$ is in agreement with frictional hydraulic theory that states that frictional effects reduce the maximum possible value of $Q_m$ (Gu & Lawrence Reference Gu and Lawrence2005). It is only for $\theta \geq 1.0^\circ$ that the additional forcing by gravity in the along-duct direction reduces the importance of frictional effects and makes the direct transition from the HGV-A approximation to the hydraulic limit ($Q_m=1/2$) possible.

Figure 8.

Mass flux per unit width $Q_m$ as a function of the parameter ${{Fr}}^*/2$ for the experiments in the mSID set-up. The colour denotes the value of the angle $\theta$. The solid black lines represent $Q_m={{Fr}}^*/2$ (the expected value in the HGV-A approximation) and $Q_m=1/2$ (the expected value in the hydraulic limit). The experimental values tend to follow the trend of the prediction based on the HGV-A approximation $Q_m\approx Fr^*/2$ until ${{Fr}}^*/2\approx 0.5$. For larger values of ${{Fr}}^*/2$, a maximum, constant value of $Q_m$ is reached for each value of $\theta$ as expected from (frictional) two-layer hydraulics. The overall maximum value is $Q_m\approx 0.5$ as predicted for the hydraulic limit.

Figure 8 is reminiscent of similar plots by Hogg et al. (Reference Hogg, Ivey and Winters2001) where they plotted $Q_m$ as a function of $({{Re}}_g A^{-1})^2$ for simulations with $Sc=1$ in a horizontal duct. In that case, they observed the change from the VAD solution to the hydraulic limit. In our case with $Sc\approx 700$, we observe a similar change but from the HGV-A approximation to the hydraulic limit. In short, the flow changes from the HGV-A approximation where the value of $K$ is constant for a given value of $\theta /\alpha$, to a hydraulically controlled flow where the value of ${{Fr}}$ is constant for a given value of $\theta /\alpha$.

In general, the flows in the ${\textsf {L}}$ and ${\textsf {H}}$ regimes are very well described by the HGV-A approximation. For the ${\textsf {L}}$ experiment, ${{Fr}}=0.40$, ${{Re}}_g=398$ and $A=30$, meaning that ${{Fr}}{{Re}}_g A^{-1} \approx 5.3$ or $Re_g\approx 13.3 A$. We have seen in figure 6 that inertia can be neglected for this experiment. Furthermore, if we consider the series of experiments with $\theta =1^\circ$ ($K\approx 62$) shown in figure 8, we see that $Q_m\approx 0.5$ for ${{Fr}}^*/2>0.5$. This suggests that the limit of validity of the HGV-A (in which inertia is neglected) is around ${{Fr}}^*\approx 1$. For this set of experiments, ${{Fr}}^*=1$ corresponds to ${{Re}}_g\approx 62 A$. These two examples show that the restrictive condition ${{Re}}_g\ll A$ derived using a scaling analysis by Lefauve & Linden (Reference Lefauve and Linden2020a) to neglect advection terms can be drastically relaxed.

4.3. Regime transitions

In this section we finally discuss how the values of ${{Fr}}^*$ are related to regime transitions. Figure 9 shows the location of the different regimes in the parameter space $(\theta,{{Re}}_g)$ for all experimental set-ups with superimposed curves of constant ${{Fr}}^*$. It can be clearly seen that these curves describe the transitions between the different regimes through all values of $\theta /\alpha$. For $\theta /\alpha \gtrsim 1$, curves given by ${{Fr}}^*=\text {const.}$ are equivalent to curves given by $\theta {{Re}}_g=\text {const.}$ for a given duct. This is in agreement with the transitions derived for forced flows by Lefauve et al. (Reference Lefauve, Partridge and Linden2019a) using an energetics budget approach and that Lefauve & Linden (Reference Lefauve and Linden2020a) derived using a two-layer frictional hydraulics model. However, the agreement of the curves of constant ${{Fr}}^*$ value, based on the HGV-A approximation, extends to lazy flows ($\theta /\alpha <1$) including $\theta /\alpha \leq 0$ where the curves given by $\theta {{Re}}_g=\text {const.}$ do not make physical sense.

Figure 9.

Location of the different regimes in the (${{Re}}_g$,$\theta$) plane for the four different set-ups: mSID, LSID, HSID and tSID. The different symbols represent the different regimes, i.e. laminar (${\textsf {L}}$), Holmboe waves (${\textsf {H}}$), other waves (${\textsf {W}}$), intermittently turbulent (${\textsf {I}}$) and turbulent (${\textsf {T}}$). The solid lines represent curves of constant ${{Fr}}^*$ with the value indicated along the line. The dashed line represents $\theta /\alpha =1$. The dotted lines are examples of the transition between regimes for forced flows given by $\theta {{Re}}_g=\text {const.}$ as proposed by Lefauve et al. (Reference Lefauve, Partridge and Linden2019a).

In the mSID set-up, the transitions between successive regimes are given by ${{Fr}}^*\approx 0.5,$ 1 and 2. The detailed analysis of the three experiments in the previous section (shown as red open squares in figure 9a) and of the mass flux as a function of ${{Fr}}^*$ coincides with this picture. We observe the emergence of turbulence at ${{Fr}}^*\approx 1$, where we find the limit of validity of the HGV-A approximation. This transition further represents a transition to hydraulically controlled flows. It is interesting now to note that experiments in the ${\textsf {I}}$ regime (i.e. for $1 \lesssim {{Fr}}^* \lesssim 2$) and $\theta \geq 2$ have values of $Q_m$ larger than in the hydraulic limit where $Q_m\approx 1/2$. For these experiments, $Q_m$ seems to overshoot following the trend of the HGV-A approximation ($Q_m={{Fr}}^*/2$). This suggests that the flow in the ${\textsf {I}}$ regime has mixed characteristics of the HGV-A approximation and the fully turbulent, hydraulically controlled flow. Finally, the characteristics of the HGV-A approximation are negligible when ${{Fr}}^*\gtrsim 2$ and the flow is in the ${\textsf {T}}$ regime.

Surprisingly, although curves of constant ${{Fr}}^*$ match the transitions for all ducts, the particular values of ${{Fr}}^*$ associated with the transitions differ. For the LSID set-up, the values of ${{Fr}}^*$ marking these transitions are about a factor two larger, for the HSID about a factor three larger, and for the tSID about a factor five larger. In spite of the differences in these values, there is a remarkable similarity: there is always a factor two between the values of ${{Fr}}^*$ describing successive regime transitions.

The difference in the transition values for the tSID set-up with respect to all the other set-ups could have been expected since there is one non-dimensional parameter ($B$) that is different. In fact, the tSID set-up is much thinner than the other set-ups and, hence, 3-D effects are more important. In particular, there is additional friction from the side walls, which results in a smaller value of ${{Fr}}$ for a given forcing magnitude. This, in turn, means that the successive transitions are delayed. Further investigation of the 3-D effects would require extension of the 2-D solution presented in § 3.3 to three dimensions, which we believe is possible, but outside the scope of the current paper.

On the other hand, the difference between the mSID and LSID (both with $A=30$), and the HSID (with $A=15$) set-ups needs further explanation since the HGV-A approximation accounts for having different $A$ values. It could be that $A=15$ is not long enough for the limit of long ducts to hold. However, we do not have enough information to test this hypothesis. Still, the most surprising difference is between the LSID and mSID set-ups, as already mentioned by Lefauve & Linden (Reference Lefauve and Linden2020a), since both set-ups have the same values for $A$, $B$ and ${{Sc}}$. If the values of these non-dimensional parameters are the same, the fact that two experiments with the same values of ${{Re}}_g$ and $\theta$ are in different regimes would go against the principle of dynamic similarity. Lefauve & Linden (Reference Lefauve and Linden2020a) already mentioned that there might be a missing length scale (and its accompanying non-dimensional parameter), but they could not think of a relevant one. If one is to consider just the duct itself, there are no other non-dimensional parameters that are needed to describe the problem. For example, the non-dimensional parameters used so far are the only parameters needed to set-up a numerical simulation such as those by Kaptein et al. (Reference Kaptein, van de Wal, Kamp, Armenio, Clercx and Duran-Matute2020). This would suggest that there are certain characteristics of the set-up outside the duct that can delay or stimulate flow transitions. Unfortunately, there is much less data available for the LSID and HSID set-ups when compared with the mSID set-up. For example, there are no 3-D velocity and density fields, and there are fewer experiments for which the value of $Q_m$ is available. However, we saw that, for the mSID set-up, the different data give a self-consistent picture. Ultimately, new independent experimental or numerical efforts will be needed to clarify the differences between set-ups.

We believe that the key to understand further the regime transitions, and possibly, the differences between the mSID and LSID set-ups lies at the edges of the duct. There are three main reasons for this: (i) the flow is largely determined by flow details at the edges of the duct (e.g. the analytical solutions obtained using either two-layer hydraulics or the HGV-A approximation are highly dependent on the boundary conditions there); (ii) the HGV-A approximation formally breaks down immediately outside the duct; (iii) instabilities first emerge at the edges of the duct where the nonlinear terms are larger since both $u$ and $\partial u /\partial x$ grow as $|x|\to 1$ (this can be deduced from conservation of mass and the shape of the interface). Hence, further knowledge of the flow and the experimental set-ups at the edges of the duct is needed to clarify these points. Note, for example, that Meyer & Linden (Reference Meyer and Linden2014) used two different shapes of edges in the duct, and any of these shapes are quite different from the solid long vertical walls in the simulations by Kaptein et al. (Reference Kaptein, van de Wal, Kamp, Armenio, Clercx and Duran-Matute2020). Furthermore, more information about the flow at the edges of the duct can help us understand how it transitions from the HGV-A approximation to being hydraulically controlled.

In spite of the differences between the set-ups, figure 9 suggests that the parameter ${{Fr}}^*$ does quantify the forcing for a given set-up. Recall that this parameter groups together the control parameters: ${{Re}}_g$, $A$ and $\theta$ with the latter incorporated in the value of $K$. It was already known that ${{Re}}_g A^{-1}$ is the control parameter for horizontal ducts (Hogg et al. Reference Hogg, Ivey and Winters2001). In this way, ${{Fr}}^*$ is an equivalent number that is only rescaled by dividing by $K$ to account for the relative inclination of the duct with respect to its internal angle (i.e. depending on $\theta /\alpha$). It would seem that $K$ is an intrinsic parameter of SID set-ups. A different way to analyse the problem is in terms of a generalized gravitational Reynolds number for inclined ducts. If the ${{Re}}_g$ quantifies the forcing due to the density differences, a generalized gravitational Reynolds number

(4.3)\begin{equation} {{Re}}_{g,\theta}= \frac{K(\theta=0)}{K}{{Re}}_g\approx \frac{131}{K} {{Re}}_g \end{equation}

can be defined, in a similar way as done by Kaptein (Reference Kaptein2021). In this way, ${{Re}}_{g,\theta }$ incorporates both aspects of the forcing: the density difference between the tanks and the tilting of the duct.

We have shown that the transition from laminar to turbulent in a (long) SID experiment with ${{Sc}}\gg 1$ depends solely on ${{Fr}}^*$. For weak forcing (small ${{Fr}}^*$ values), the flow is laminar and can be described by HGV-A approximation. As the forcing (the value of ${{Fr}}^*$) is increased, the flow develops nonlinear and time-dependent characteristics (waves, thinner boundary layers, turbulence). Once the flow becomes turbulent, the HGV-A approximation becomes invalid. The nonlinear characteristics result in additional vertical viscous momentum diffusion that, for ${{Fr}}^*>1$, compensates the increase in forcing. In this way, the value of ${{Fr}}$ can remain constant in agreement with (frictional) two-layer hydraulic theory (Gu & Lawrence Reference Gu and Lawrence2005).