Inhomogeneous Couette flows for a two-layer fluid

Full Text

Introduction

The Couette flow is one of the first examples of an exact solution to the Navier–Stokes equations [1–3]. It describes the isobaric flow of a viscous incompressible fluid, which is induced due to the motion of one or both boundaries of an infinite horizontal fluid layer. Recall that the steady Couette flow is described by a linear velocity profile, and this has predetermined the popularity of this solution in the theory of hydrodynamic stability for studying secondary flows generated by different disturbance classes. 
There is an interesting observation for the unsteady Couette flow. It is described by the simplest linear parabolic equation having a general solution to an extensive class of functions [2, 3]. Thus, the first and second Stokes problems are described by the non-stationary Couette profile, and they are its particular case [2, 3].

Besides studying the hydrodynamic stability of the Couette flow for a viscous incompressible fluid, there are various modifications of the exact Couette solution for regions without plane symmetry. Note that there exists the well-known exact Taylor–Couette solution on the fluid flow in the gap between coaxial cylinders [2, 4–8], as well as the solution describing isobaric fluid flow on a sphere [9, 10]. The Couette flow is fundamental in the study of fluids with non-Newtonian properties [11, 12]. The three-dimensional Couette flow, which is potential and non-isobaric, has been recently exemplified [13].

It is difficult to study Couette flows, different from unidirectional ones, since the reduced Navier–Stokes system becomes overdetermined [14–17]. The overdetermined Navier–Stokes equation system describing two-dimensional flows of viscous fluids began to be studied in [17], a complete list of exact solutions for two-dimensional hydrodynamics being given in [18]. An example of a nontrivial exact solution to the Navier–Stokes equation system for incompressible fluids with nonstationary and steady two-dimensional velocity fields depending on three coordinates was found in [15, 16]. The first exact solution describing the non-uniform Couette flow was constructed in the class of solutions for velocities linearly dependent on two coordinates (the Lin–Sidorov–Aristov family) [19–21]. The studies along this line were continued and summarized in [16].

This paper studies a boundary value problem for 2.5D Navier–Stokes equations, which describes steady flows of a two-layer fluid. The study is based on the exact solution of the Navier–Stokes equations for incompressible fluids [15, 16], which was constructed by functional variable separation. It was shown in [22] that the exact solutions found in [15, 16]
and describing non-uniform Couette flows can be used to describe isobaric multilayer fluids. It was reported in [15, 16] that exact solutions with a velocity field linear in coordinates describe equatorial countercurrents in the World Ocean [2, 15, 21]. Multilayer fluids are often used to model large oceanic flows [2, 21]. In [22] it was found useful to extend the study of boundary value problems of steady flows from single-layer streams [15, 16] to multilayer fluids. This should be done in order to have a store of exact solutions for studying the hydrodynamic stability of flows, for comparing model representations with natural observations, and most importantly to understand the applicability of the substitution of continuous density stratification by a discrete saltus function.

1. Problem Statement

Consider a steady flow of a viscous two-layer fluid in an extended horizontal layer. We denote the lower layer by the subscript "1" and the upper one by "2". Each of the two layers of the two-layer fluid can have its own thickness ($h_1$ and $h_2$, respectively, see Fig. 1), density ($\rho_1$ and $\rho_2$), and viscosity ($\eta_1$ and $\eta_2$). Note that the heavier phase is located below in density-stratified fluids, i.e., $\rho_1 > \rho_2$.

Figure. 1. The scheme of a two-layer fluid

It is assumed here that the flow takes place at a constant temperature and in the absence of external forces, except for gravity. The sum of pressure (ratioed to density) and the gravity potential is constant through the flow; therefore, the gradient of this sum, entering the Navier–Stokes equations, is zero. The steady flow of fluids of this type is described by two systems of nonlinear equations (each layer is described by its own system) consisting of the vector Navier–Stokes equation and the scalar incompressibility equation [22–24]:
\[ \begin{equation}
 {\rho _1} ( {\mathbf{V}^{(1)}}, \nabla ){\mathbf{V}^{(1)}}={\eta_1}\Delta {\mathbf{V}^{(1)}}, \quad \nabla\cdot\mathbf{V}^{(1)}=0;
\end{equation} \tag{1} \]
\[ \begin{equation}
 {\rho _2} ( {\mathbf{V}^{(2)}}, \nabla ){\mathbf{V}^{(2)}}={\eta_2}\Delta {\mathbf{V}^{(2)}}, \quad \nabla\cdot\mathbf{V}^{(2)}=0.
\end{equation} \tag{2} \]

Systems (1), (2) have the following notations: $\mathbf{V}^{(1)} = (V_x^{(1)}, V_y^{(1)}, V_z^{(1)})$ and $\mathbf{V}^{(2)} = (V_x^{(2)}, V_y^{(2)}, V_z^{(2)})$ are the vector velocity fields for the lower and upper layers, respectively; $\nabla = (\partial/\partial x, \partial/\partial y, \partial/\partial z )$ is the Hamilton operator; $\Delta = (\partial^2/\partial x^2 + \partial^2/\partial y^2 + \partial^2/\partial z^2 )$ is the Laplace operator; $(\mathbf{V},\nabla ) = (V_x \partial/\partial x + V_y \partial/\partial y + V_z \partial/\partial z )$ is a convective derivative. Note that both system (1) and system (2) are overdetermined, i.e., the required velocities $V_x^{(i)}$, $V_y^{(i)}$ for each layer must satisfy three scalar equations (the third Navier–Stokes equation is fulfilled identically since a flow with zero vertical velocity is considered). Note that the equations $\nabla\cdot\mathbf{V}^{(1)}=0$ and $\nabla\cdot\mathbf{V}^{(2)}=0$ (velocity divergence is zero) is termed in two ways in the scientific literature. They are termed the continuity equation in the physical literature and the incompressibility equation in the hydrodynamic literature [25, 26].

In what follows, the solution of systems (1), (2) is sought in the form of functions linearly dependent on one of the horizontal coordinates [22–24]:
\[ \begin{equation}
 {V_x}^{(1)}={u^{(1)}} ( z )+{u^{(1)}} ( z )y, \quad {V_y}^{(1)}={V^{(1)}} ( z );
\end{equation} \tag{3} \]
\[ \begin{equation}
 {V_x}^{(2)}={u^{(2)}} ( z )+{u^{(2)}} ( z )y, \quad {V_y}^{(2)}={V^{(2)}} ( z ).
\end{equation} \tag{4} \]

For the class (3), (4) the incompressibility equation in systems (1), (2) is satisfied identically. This circumstance dismantles the problem of overdetermination of systems (1) and (2), i.e., each system is now reduced to finding two projections of the velocity vector from two ordinary differential equations (projections of the Navier–Stokes equation onto the axes $Ox$ and $Oy$):
\[ \begin{equation}
 {\rho _1}{V^{(1)}}{u^{(1)}}={\eta_1} ( {u^{(1)}}''+ {u^{(1)}}''y ), \quad \eta _{1} {V^{(1)}}''=0;
\end{equation} \tag{5} \]
\[ \begin{equation}
 {\rho _1}{V^{(2)}}{u^{(2)}}={\eta_2} ( {u^{(2)}}''+{u^{(2)}}''y ), \quad \eta _{2} {V^{(2)}}''=0.
\end{equation} \tag{6} \]
Hereinafter, the double prime marks derivation with respect to the vertical coordinate $z$. In view of the independence of the spatial coordinates $x$ and $y$ of the selected Cartesian system, equations (5) and (6) can be represented as
\[ \begin{equation}
 {u^{(1)}}''=0, \quad {V^{(1)}}''=0, \quad {u^{(1)}}''=\frac{{\rho _1}}{{\eta_1}}{V^{(1)}}{u^{(1)}};
\end{equation} \tag{7} \]
\[ \begin{equation}
 {u^{(2)}}''=0, \quad {V^{(2)}}''=0, \quad {u^{(2)}}''=\frac{{\rho _2}}{{\eta_2}}{V^{(2)}}{u^{(2)}}.
\end{equation} \tag{8} \]
The first two equations in both system (7) and system (8) are isolated, and the solution of the third equations in these systems is the last to be found. Double integration of systems (7) and (8) results in their general solution
\[ \begin{equation*}
{u^{(1)}}={{c}_{2}}^{(1)}z+{{c}_{1}}^{(1)}, \quad {V^{(1)}}={{\alpha }_{2}}^{(1)}z+{{\alpha }_{1}}^{(1)},
\end{equation*} \]
\[ \begin{equation}
{u^{(1)}}=\frac{{\rho _1}}{12{\eta_1}} 
\bigl( {{c}_{2}}^{(1)}{{\alpha }_{2}}^{(1)}{{z}^{4}}+2 ( {{c}_{2}}^{(1)}{{\alpha }_{1}}^{(1)}+{{c}_{1}}^{(1)}{{\alpha }_{2}}^{(1)} ){Z^3}+6{{c}_{1}}^{(1)}{{\alpha }_{1}}^{(1)}{{z}^{2}} \bigr)+{{\beta }_{2}}^{(1)}z+{{\beta }_{1}}^{(1)};
\end{equation} \tag{9} \]
\[ \begin{equation*}
{u}^{(2)}={{c}_{2}}^{(2)}z+{{c}_{1}}^{(2)}, \quad {V^{(2)}}={{\alpha }_{2}}^{(2)}z+{{\alpha }_{1}}^{(2)},
\end{equation*} \]
\[ \begin{equation}
{u^{(2)}}=\frac{{\rho _2}}{12{\eta_2}}
 \bigl( {{c}_{2}}^{(2)}{{\alpha }_{2}}^{(2)}{{z}^{4}}+2 
 ( {{c}_{2}}^{(2)}{{\alpha }_{1}}^{(2)}+{{c}_{1}}^{(2)}{{\alpha }_{2}}^{(2)} ){Z^3}+6{{c}_{1}}^{(2)}{{\alpha }_{1}}^{(2)}{{z}^{2}} \bigr)+{{\beta }_{2}}^{(2)}z+{{\beta }_{1}}^{(2)}.
\end{equation} \tag{10} \]
The solutions represented by equations (9) and (10) are polynomial, the highest degree of these polynomials corresponds to expressions for velocities $U^{(1)}$, $U^{(2)}$, and this is attributable to the sequence of integration of the equations in systems (7) and (8). The constants $c_i^{(1)}$, $c_i^{(2)}$, $\alpha_i^{(1)}$, $\alpha_i^{(2)}$, $\beta_i^{(1)}$ and $\beta_i^{(2)}$ $(i = 1, 2 )$ in the exponential solutions (9) and (10) must be found from the boundary conditions; therefore, it is necessary to formulate twelve conditions for the determination of these values.

2. Boundary Conditions

Since the vertical fluid velocity is assumed to be zero, the fluids of the different layers do not intermix in the shear flow under study. In other words, the interlayer boundary (the boundary $z = h_1$) is here considered to be rigid (Fig. 1). For convenience, in what follows, $h = h_1 + h_2$.

Assume that a no-slip condition is set at the lower boundary $z = 0$ [15]:
\[ \begin{equation*}
V_x^{(1)} (0 )=0, \quad V_y^{(1)} (0 )=0.
\end{equation*} \]

Taking into account the representation (3), (4), we arrive at three conditions
\[ \begin{equation} 
 U^{(1)} (0 )=0, \quad u^{(1)} (0 )=0, \quad V^{(1)} (0 )=0.
\end{equation} \tag{11} \]

A velocity field (wind effect) is set at the upper boundary $z = h$ [15]:
\[ \begin{equation*}
V_x^{(2)}(h)=W\cos\varphi+\Omega y, \quad V_y^{(2)}(h)=W\sin\varphi.
\end{equation*} \]
Here, $\Omega$ is the horizontal gradient of the velocity $V_x^{(2)}$ (spatial acceleration) at the upper boundary; $W$ is the absolute value of the uniform velocity component at the mobile boundary $z = h$ of a two-layer fluid; $\varphi$ is the angle between the uniform velocity component $V_x^{(2)}$ and the axis $Ox$. Taking into account representation (4), we obtain three more conditions
\[ \begin{equation}
 U^{(2)}(h)=W\cos\varphi, \quad u^{(2)}(h)=\Omega, \quad V^{(2)}(h)=W\sin\varphi.
\end{equation} \tag{12} \]

Besides, it is required that two additional conditions be met at the interlayer boundary $z = h_1$ (conditions for "sewing together" the solutions for the two layers). They are

• the solution continuity condition for velocities 
  \[ \begin{equation*}
  V_x^{(1)}(h_1)=V_x^{(2)}(h_1), \quad 
  V_y^{(1)}(h_1)=V_y^{(2)}(h_1),
  \end{equation*} \]
  and taking into account the structure of classes (3) and (4), we have three equalities
  \[ \begin{equation}
  U^{(1)}(h_1)=U^{(2)}(h_1), \quad u^{(1)}(h_1)=u^{(2)}(h_1), \quad
  V^{(1)}(h_1)=V^{(2)}(h_1);
  \end{equation} \tag{13} \]
• the solution continuity condition for tangential stresses
  \[ \begin{equation*}
  \tau_{xz}^{(1)}(h_1)=\tau_{xz}^{(2)}(h_1), \quad \tau_{yz}^{(1)}(h_1)=\tau_{yz}^{(2)}(h_1).
  \end{equation*} \]

As distinct from expressions (13) resulting (in view of the independence of the spatial coordinates) directly from the condition of equality of the velocities at the interlayer boundary, the case with the continuity condition for tangential stresses is not as apparent, the continuity condition for the velocity field proves to be insufficient. By the Newton law, the relation of the stress tensor components to velocities is known,
\[ \begin{equation*}  {\tau ^{(i)}}=\left(   \begin{matrix} -{p^{(i)}}+2{\eta_i}\frac{\partial {V_x}^{(i)}}{\partial x}  &  {\eta_i}\bigl( \frac{\partial {V_x}^{(i)}}{\partial y}+\frac{\partial {V_y}^{(i)}}{\partial x} \bigr)  &  {\eta_i}\bigl( \frac{\partial {V_x}^{(i)}}{\partial z}+\frac{\partial {V_z}^{(i)}}{\partial x} \bigr)  \\ {\eta_i}\bigl( \frac{\partial {V_x}^{(i)}}{\partial y}+\frac{\partial {V_y}^{(i)}}{\partial x} \bigr)  &  -{p^{(i)}}+2{\eta_i}\frac{\partial {V_y}^{(i)}}{\partial y}  &  {\eta_i}\bigl( \frac{\partial {V_y}^{(i)}}{\partial z}+\frac{\partial {V_z}^{(i)}}{\partial y} \bigr)  \\ {\eta_i}\bigl( \frac{\partial {V_x}^{(i)}}{\partial z}+\frac{\partial {V_z}^{(i)}}{\partial y} \bigr)  &  {\eta_i}\bigl( \frac{\partial {V_z}^{(i)}}{\partial y}+\frac{\partial {V_y}^{(i)}}{\partial z} \bigr)  &  -{p^{(i)}}+2{\eta_i}\frac{\partial {V_z}^{(i)}}{\partial z} \\  \end{matrix}   \right)  \overset{(3), (4)}{\mathop{=}} \left(   \begin{matrix} -{p^{(i)}} & {\eta_i}{u^{(i)}} & {\eta_i}\frac{\partial {V_x}^{(i)}}{\partial z}  \\ {\eta_i}{u^{(i)}} & -{p^{(i)}} & {\eta_i}\frac{\partial {V_y}^{(i)}}{\partial z}  \\ {\eta_i}\frac{\partial {V_x}^{(i)}}{\partial z}  &  {\eta_i}\frac{\partial {V_y}^{(i)}}{\partial z}  &  -{p^{(i)}} \\ \end{matrix} \right), \end{equation*} \]
where $p^{(i)}$ is hydrostatic pressure (the pressure of a moving fluid column, depending not only on the transverse coordinate, i.e. changing only with depth, similarly to the main equation of hydrostatics) in the layer under study. This representation of the stress tensor can easily yield the following equivalent of the continuity condition for tangential stresses:
\[ \begin{equation*} 
\eta_1\frac{du^{(1)}}{dz} \Bigr|_{z=h_1}=
\eta_2\frac{du^{(2)}}{dz} \Bigr|_{z=h_1}, 
\quad
\eta_1\frac{du^{(1)}}{dz} \Bigr|_{z=h_1}=
\eta_2\frac{du^{(2)}}{dz} \Bigr|_{z=h_1}, 
\end{equation*} \]
\[ \begin{equation}
\eta_1\frac{dV^{(1)}}{dz} \Bigr|_{z=h_1}=
\eta_2\frac{dV^{(2)}}{dz} \Bigr|_{z=h_1}.
\end{equation} \tag{14} \]
Thus the boundary value problem (9)–(14) becomes closed.

3. An Exact Solution to the Boundary Value Problem

Conditions (11)–(14) allow us to find a particular solution to systems (7), (8) for each layer, which would meet the selected boundary conditions:
\[ \begin{equation}
{u^{(1)}}(Z)=\frac{h{\eta_2}\Omega }{{h_2}{\eta_1}+{h_1}{\eta_2}}Z,
\end{equation} \tag{15} \]
\[ \begin{equation}
{V^{(1)}}(Z)=\frac{h{\eta_2}W\sin\varphi }{{h_2}{\eta_1}+{h_1}{\eta_2}}Z,
\end{equation} \tag{16} \]
\[ \begin{multline}
{u^{(1)}}(Z)=\frac{hWZ}{12{\eta_1}{{ ( {h_2}{\eta_1}+{h_1}{\eta_2} )}^{3}}}
\bigl[ 
\Omega {\rho _1}{\eta_2^2} ( {h_2}{\eta_1}+{h_1}{\eta_2} ) {h^3}\sin\varphi {Z^3}+12{\eta_1}{\eta_2}{{ ( {h_2}{\eta_1}+{h_1}{\eta_2} )}^2}\cos\varphi- \\ 
-\Omega \sin\varphi ( 4{h_1^3}{h_2}{\eta_1}{\eta_2^2}{\rho _1}+{h_1^4}{\eta_2^3}{\rho _1}+{h_2^4}{\eta_1^3}{\rho _2}+ 4{h_1}{h_2^3}{\eta_2^2}{\eta_2}{\rho _2}+6{h_1^2}{h_2^2}{\eta_1}{\eta_2^2}{\rho _2} ) \bigr],
\end{multline} \tag{17} \]
\[ \begin{equation}
 {u^{(2)}}(Z)=\frac{\Omega }{{h_2}{\eta_1}+{h_1}{\eta_2}}\bigl( h{\eta_1}Z+{h_1} ( {\eta_2}-{\eta_1} ) \bigr),
\end{equation} \tag{18} \]
\[ \begin{equation}
{V^{(2)}}(Z)=\frac{W\sin\varphi }{{h_2}{\eta_1}+{h_1}{\eta_2}}\bigl( h{\eta_1}Z+{h_1} ( {\eta_2}-{\eta_1} ) \bigr),
\end{equation} \tag{19} \]
\[ \begin{multline}
{u^{(2)}}(Z)=\frac{W}{12{\eta_2}{{( {h_2}{\eta_1}+{h_1}{\eta_2} )}^{3}}}
\bigl\{ {h^4}{\eta_2^2} ( {h_2}{\eta_1}+{h_1}{\eta_2} ){\rho _2}\Omega \sin\varphi {Z^4}-4{h^3}{h_1}{\eta_1} ( {\eta_1}-{\eta_2} ) ( {h_2}{\eta_1}+{h_1}{\eta_2} ){\rho _2}\Omega \sin\varphi Z^3+ \\
+6{h^2}{h_1^2}{{ ( {\eta_1}-{\eta_2} )}^{2}} ( {h_2}{\eta_1}+{h_1}{\eta_2} ){\rho _2}\Omega \sin\varphi Z^2 -hZ\bigl[ \bigl( -3{h_1^4}{\eta_2^3}{\rho _1}+ \bigl( {h_2} ( 4{h_1^3}+{h_2^3} ){\eta_1^3}+4{h_1} ( {h_1^3}-3{h_1^2}{h_2}+{h_2^3} ){\eta_2^2}{\eta_2}+ \\
+6{h_1^2} ( -2{h_1^2}+2{h_1}{h_2}+{h_2^2} ){\eta_1}{\eta_2^2}+12{h_1^4}{\eta_2^3} \bigr){\rho _2} \bigr)
\Omega \sin\varphi -12{\eta_1}{\eta_2}{( {h_2}{\eta_1}+{h_1}{\eta_2} )}^{2}\cos\varphi \bigr]+ \\
-{h_1} ( {h_1}+{h_2} )W\sin\varphi \bigl( 3{h_1^3}{\eta_2^3}{\rho _1}-{h_1^2}{h_2}{\eta_1^3}{\rho _2}+{h_1}{h_2^2}{\eta_1^3}{\rho _2}-{h_2^3}{\eta_1^3}{\rho _2}-{h_1^3}{\eta_2^2}{\eta_2}{\rho _2}+ 5{h_1}{h_2}{\eta_2^2}{\eta_2}{\rho _2} ( {h_1}-{h_2} )+{h_2^3}{\eta_2^2}{\eta_2}{\rho _2}+4{h_1^3}{\eta_1}{\eta_2^2}{\rho _2}- \\
-10{h_1^2}{h_2}{\eta_1}{\eta_2^2}{\rho _2}+4{h_1}{h_2^2}{\eta_1}{\eta_2^2}{\rho _2}+6{h_1^2}{\eta_2^3}{\rho _2}
( {h_2}-{h_1} ) \bigr)-12{h_1} ( {\eta_1}-{\eta_2} ){\eta_2}{{( {h_2}{\eta_1}+{h_1}{\eta_2} )}^{2}}\cos\varphi.
\end{multline} \tag{20} \]
Solutions (15)–(20) involve the substitution \(Z = z/h \in [0,~1]\). The motion of the lower layer \((Z \in [0; h_1/h ])\) is described by expressions (15)–(17). The values of $Z \in [h_1/h; 1 ]$ correspond to the other (upper) layer, where the flow velocity is determined by solutions (18)–(20).

It is of interest that, despite seemingly simple structure of boundary conditions (11), (12) and velocity field representation (11), (12), the obtained solution (15)–(20) depends on the parameters of the boundary problem under study and the physical characteristics of the fluid in an extremely non-trivial way. Therefore, in a general form, it is impossible to make any conclusions about the effect of a specific parameter on the structure of the final solution. First of all, this concerns the uniform components $U^{(1)}$ and $U^{(2)}$. The only conclusion (fairly obvious) that can be made concerns the effect of the values of $W$ and $\Omega$ on the value of velocity in the $Oy$ direction (the velocities $V^{(1)}$ and $V^{(2)}$) and the non-uniform velocity component along the $Ox$ axis (the spatial gradients $u^{(1)}$ and $u^{(2)}$) since functions (15), (16) and (18), (19) depend in direct proportion on these parameters.

Note also that, when $W = 0$, the velocity field determined by solutions (15)–(17) and (18)–(20) assumes a trivial form; therefore, it is considered hereinafter that $W \neq 0$.

Besides, note that the functions $u^{(1)}$ and $u^{(2)}$ become zero when $\Omega = 0$, i. e., we obtain an extension of the classical Couette flow to a non-one-dimensional case. The velocities $V^{(1)}$ and $V^{(2)}$ become zero only if $\sin\varphi = 0$. In this case, there is a unidirectional flow along the $Ox$ axis with a non-uniform velocity distribution. In addition, the exact solutions for the components $u^{(1)}$ and $V^{(1)}$ are linear functions of $Z$, whose zero is only the point $Z = 0$. The stagnation (zero) points of the fluid are of particular importance due to the fact that both simply flow suppression and flow reversal (the appearance of counterflows) are possible at these points.

4. Analysis of the Exact Solution for the Lower Layer

Expression (17) for velocity $U^{(1)}$ can be written as
\[ \begin{equation}
{u^{(1)}}=\frac{hWZ}{12{\eta_1}{{( {h_2}{\eta_1}+{h_1}{\eta_2} )}^{3}}} [ a{Z^3}+b ],
\end{equation} \tag{21} \]
where the following nomenclature is introduced:
\[ \begin{equation*}
a=\Omega {\rho _1}{\eta_2^2} ( {h_2}{\eta_1}+{h_1}{\eta_2} ){h^3}\sin\varphi,
 \end{equation*} \]
\[ \begin{equation*}
b=12{\eta_1}{\eta_2}{{ ( {h_2}{\eta_1}+{h_1}{\eta_2} )}^{2}}\cos\varphi - \Omega \sin\varphi 
( 4{h_1^3}{h_2}{\eta_1}{\eta_2^2}{\rho _1}+{h_1^4}{\eta_2^3}{\rho _1}
 +{h_2^4}{\eta_1^3}{\rho _2}+4{h_1}{h_2^3}{\eta_2^2}{\eta_2}{\rho _2}+6{h_1^2}{h_2^2}{\eta_1}{\eta_2^2}{\rho _2} ).
\end{equation*} \]
The case $a = 0$ in Eq. (21) is not discussed here due to its triviality. Assume further that $a \neq 0$, then the velocity $U^{(1)}$ becomes zero inside the layer under study only if the point \(Z_0 = - \sqrt[{3}]{b/a}\) falls within this layer. The latter condition is equivalent to the fulfillment of the inequality
\[ \begin{equation*}
b (ah_1^3+bh^3 )<0.
\end{equation*} \]
The corresponding velocity profiles $U^{(1)}$ are shown in Figs. 2 and 3. These figures show the value of angular velocity, corresponding to the first Coriolis parameter, for the equatorial zone of the World Ocean (the latitude is $\pi/12$ rad).

Figure. 2. The \(U\)-velocity profiles (the dashed line \(U^{(1)}\) for the lower layer and the solid line \(U^{(2)}\) for the upper layer) for \(W=10\) m/s, \(\varphi=\pi/12\), \(h_1=0.4\) m, \(h_2=0.6\) m, \(\Omega1.4584\cdot10^{-5}\) rad/s, \(\rho_1=1100\) kg/m\(^3\), \(\rho_2=1057.6\) kg/m\(^3\), \(\eta_1=1.7\cdot10^{-5}\) Pa\(\cdot\)s, \(\eta_2=1.2\cdot10^{-5}\) Pa\(\cdot\)s

Figure. 3. The \(U\)-velocity profiles (the dashed line \(U^{(1)}\) for the lower layer and the solid line \(U^{(2)}\) for the upper layer) for \(W=10\) m/s, \(\varphi=\pi/12\), \(h_1=0.6\) m, \(h_2=0.4\) m, \(\Omega1.4584\cdot10^{-5}\) rad/s, \(\rho_1=1100\) kg/m\(^3\), \(\rho_2=733.1\) kg/m\(^3\), \(\eta_1=1.7\cdot10^{-5}\) Pa\(\cdot\)s, \(\eta_2=1.2\cdot10^{-5}\) Pa\(\cdot\)s

In view of the above notations, the tangential stress $\tau_{xz}^{(1)}$ is defined as
\[ \begin{equation*}
\tau _{xz}^{(1)}=\frac{W}{12{{ ( {h_2}{\eta_1}+{h_1}{\eta_2} )}^{3}}} [ 4a{Z^3}+b ]+\frac{{\eta_1}{\eta_2}\Omega y}{{h_2}{\eta_1}+{h_1}{\eta_2}}.
\end{equation*} \]
Note that, when $W = 0$, the tangential stress $\tau_{xz}^{(1)}$ assumes a constant value (different for each section $y$). Therefore, it is further assumed that $W \neq 0$, and hence the stress $\tau_{xz}^{(1)}$ can be represented as
\[ \begin{equation}
 \tau _{xz}^{(1)}= \frac{hW}{12{{( {h_2}{\eta_1}+{h_1}{\eta_2} )}^{3}}} [ 4a{Z^3}+b ]+\frac{h{\eta_1}{\eta_2}\Omega y}{{h_2}{\eta_1}+{h_1}{\eta_2}}=\frac{hW}{12{{( {h_2}{\eta_1}+{h_1}{\eta_2} )}^{3}}} [ 4a{Z^3}+C ],
\end{equation} \tag{22} \]
where $C$ denotes the following expression:
\[ \begin{equation*}
 C=b+\frac{12{\eta_1}{\eta_2}{{ ( {h_2}{\eta_1}+{h_1}{\eta_2} )}^{2}}\Omega y}{W}.
\end{equation*} \]

Note that the structure of expression (22) is similar to formula (21). Hence, by analogy with the analysis of the expression (21), it follows that, when $a \neq 0$, the stress $\tau_{xz}^{(1)}$ becomes zero only at one point inside the layer, namely at the point \(Z = - \sqrt[{3}]{C/(4a)}\). This value belongs to the interval of interest only if
\[ \begin{equation*}
 C (4ah_1^3+Ch^3 )<0.
\end{equation*} \]
Thus, the profile of the tangential stress represented by expression (22) is defined by a cubic parabola, which cannot have more than one point in common with the $OZ$ axis (Figs. 4 and 5). The profiles of the tangential stress $\tau_{xz}$ in Figs. 4 and 5 are constructed for those values of the boundary conditions, layer thicknesses, and physical fluid parameters for which a calculation is made and presented in Figs. 2 and 3, respectively.

Figure. 4. The profile of the stress \(\tau_{xz}\) (the dashed line \(\tau_{xz}^{(1)}\) for the lower layer and the solid line \(\tau_{xz}^{(2)}\) for the upper layer) for \(W=10\) m/s, \(\varphi=\pi/12\), \(h_1=0.4\) m, \(h_2=0.6\) m, \(\Omega1.4584\cdot10^{-5}\) rad/s, \(\rho_1=1100\) kg/m\(^3\), \(\rho_2=1057.6\) kg/m\(^3\), \(\eta_1=1.7\cdot10^{-5}\) Pa\(\cdot\)s, \(\eta_2=1.2\cdot10^{-5}\) Pa\(\cdot\)s

Figure. 5. The profile of the stress \(\tau_{xz}\) (the dashed line \(\tau_{xz}^{(1)}\) for the lower layer and the solid line \(\tau_{xz}^{(2)}\) for the upper layer) for \(W=10\) m/s, \(\varphi=\pi/12\), \(h_1=0.6\) m, \(h_2=0.4\) m, \(\Omega1.4584\cdot10^{-5}\) rad/s, \(\rho_1=900\) kg/m\(^3\), \(\rho_2=733.1\) kg/m\(^3\), \(\eta_1=1.7\cdot10^{-5}\) Pa\(\cdot\)s, \(\eta_2=1.2\cdot10^{-5}\) Pa\(\cdot\)s

Note that, in view of the structure of solutions (16) and (19), the tangential stress $\tau_{yz}$ assumes a constant value determined by the problem parameters:
\[ \begin{equation*}
\tau _{yz}^{(1)}={\eta_1}\frac{dV^{(1)}}{dz}=\frac{{\eta_1}}{h} \frac{dV^{(1)}}{dZ}=\frac{{\eta_1}{\eta_2}W \sin\varphi }{{h_2}{\eta_1}+{h_1}{\eta_2}}={{\tau }_{yz}}^{(2)}.
\end{equation*} \]

5. Analysis of the Exact Solution for the Upper Layer

The expression (20) for the velocity $U^{(2)}$ can be represented as
\[ \begin{equation}
 {u^{(2)}}=\frac{W}{12{\eta_2}{ ( {h_2}{\eta_1}+{h_1}{\eta_2} )^{3}}} [ k{{Z}^{4}}+n{Z^3}+m{{Z}^{2}}+qZ+r ].
\end{equation} \tag{23} \]
The coefficients $k$, $n$, $m$, $q$ and $r$ can be easily written from exact solution (20).

Let us now study the properties of the tangential stress defined by the velocity $V_x^{(2)}$. Using the above notations, we obtain
\[ \begin{equation*}
\tau _{xz}^{(2)}=\frac{{\eta_2}}{h} \Bigl\{ \frac{W}{12{\eta_2}{{( {h_2}{\eta_1}+{h_1}{\eta_2} )}^{3}}}
[ 4k{Z^3}+3n{{Z}^{2}}+2mZ+q ]+\frac{h{\eta_1}\Omega y}{{h_2}{\eta_1}+{h_1}{\eta_2}} \Bigr\}.
\end{equation*} \]
Note that, when $W = 0$, the stress $\tau_{xz}^{(2)}$ assumes a constant (through the selected section) value, hence the stratification of this field does not occur. In the assumption that $W \neq 0$, the stress $\tau_{xz}^{(2)}$ can be represented as
\[ \begin{equation*}
\tau _{xz}^{(2)}=\frac{{\eta_2}}{h}\frac{W}{12{\eta_2}{{ ( {h_2}{\eta_1}+{h_1}{\eta_2} )}^{3}}} \Bigl\{ 4k{Z^3}+3n{{Z}^{2}}+2mZ+q + \frac{12h{\eta_1}{\eta_2}{{ ( {h_2}{\eta_1}+{h_1}{\eta_2} )}^{2}}\Omega y}{W} \Bigr\},
\end{equation*} \]
and this illustrates stratification on the change in the sign of the tangential stress relative to the $Z$-coordinate. The examples of the tangential stress $\tau_{xz}^{(2)}$ corresponding to this solution are given in Figs. 4 and 5.

Conclusion

The problem describing the isothermal flow of a viscous incompressible two-layer fluid in a horizontal layer has been studied. The properties of the layers differ in thickness, density and/or viscosity. An exact solution has been obtained to describe velocities in each of the layers for a set of boundary conditions describing no-slip of the fluid at the lower boundary of the fluid flow region under study and the non-uniform effect of wind at the upper boundary of this region. At the boundary between the layers of the two-layer fluid, it was required that the velocities and stresses be equal. The analysis of the obtained solution for the layers has shown that this exact solution is able to describe the appearance of reverse flow zones and stratification of the tangential stress field.

Competing interests. Authors declare that they have no competing interests.
Authors’ contributions and responsibilities. Each author has participated in the development of the concept of the article and in the writing of the manuscript. The authors are absolutely responsible for submitting the final manuscript in print. Each author has approved the final version of the manuscript.
Funding. The study received no funding.

About the authors

Natalya V. Burmasheva

Institute of Engineering Science, Urals Branch, Russian Academy of Sciences; Ural Federal University named after the First President of Russia B. N. Yeltsin

Author for correspondence.
Email: nat_burm@mail.ru
ORCID iD: 0000-0003-4711-1894
SPIN-code: 7927-5530
Scopus Author ID: 57193346922
ResearcherId: E-3908-2016
http://www.mathnet.ru/person52636

Senior Researcher, Sect. of Nonlinear Vortex Hydrodynamics; Cand. Tech. Sci., Associate Professor, Dept. of Information Technology and Automation

Russian Federation, 620049, Ekaterinburg, Komsomolskaya st., 34; 620002, Ekaterinburg, Mira st., 19

Ekaterina A. Larina

Ural Federal University named after the first President of Russia B. N. Yeltsin; Institute of Engineering Science, Ural Branch of RAS

Email: larinakatia@yandex.ru
ORCID iD: 0009-0000-7883-0803

Assistant, Dept. of Information Technology and Automation; Engineer, Sect. of Nonlinear Vortex Hydrodynamics

Russian Federation, 620002, Ekaterinburg, Mira st., 19; 620049, Ekaterinburg, Komsomolskaya st., 34

Evgeniy Yu. Prosviryakov

Institute of Engineering Science, Urals Branch, Russian Academy of Sciences; Ural Federal University named after the First President of Russia B. N. Yeltsin; Udmurt Federal Research Center, Ural Branch of RAS

Email: evgen_pros@mail.ru
ORCID iD: 0000-0002-2349-7801
SPIN-code: 3880-5690
Scopus Author ID: 57189461740
ResearcherId: E-6254-2016
http://www.mathnet.ru/person41426

Head of Sector, Sect. of Nonlinear Vortex Hydrodynamics; Dr. Phys. & Math. Sci., Professor, Dept. of Information Technology and Automation; Leading Researcher, Lab. of Physical and Chemical Mechanics

Russian Federation, 620049, Ekaterinburg, Komsomolskaya st., 34; 620002, Ekaterinburg, Mira st., 19; 426067, Izhevsk, T. Baramzina st., 34

References

Supplementary files