2. Mathematical modelling

We consider a liquid film flowing down a plane making an angle $\beta$ with the horizontal. The fluid properties, i.e. kinematic viscosity $\nu =\mu /\rho$, density $\rho$ and surface tension $\sigma$, are constant. At the free surface, the gas is considered to be passive, with a constant pressure. We assume a slow evolution of the liquid film in space and time and introduce a film parameter $\varepsilon$ as a formal parameter, which is inserted along with each derivative $\partial _x$ and $\partial _t$ in space and time. This allows us to apply Prandtl's simplification of the cross-stream momentum balance, which gives

(2.1a)\begin{gather} \varepsilon\delta (u_t + u u_x + v u_y) = b(h) + u_{yy} + \varepsilon^2 \eta (2u_{xx} + [ u_x|_{y=h} ]_x), \end{gather}

(2.1b)\begin{gather}u_x + v_y = 0, \quad u|_{y=0}= v|_{y=0}, \end{gather}

(2.1c)\begin{gather}u_y|_{y=h} = \varepsilon^2\eta (4 h_x u_x|_{y=h} - v_x|_{y=h}),\quad v|_{y=h} = h_t + u|_{y=h} h_x. \end{gather}

Here $b(h) = 1 - \varepsilon \zeta h_x + \varepsilon ^3 h_{xxx}$ combines the streamwise gravity acceleration (body force) and a pressure gradient with a hydrostatic part and a contribution of surface tension.

System (2.1) is obtained after elimination of the pressure field and is consistent up to $O(\varepsilon ^2)$ (Ruyer-Quil & Manneville Reference Ruyer-Quil and Manneville2000). We chose to write (2.1) using Shkadov's notation, for which the time $t$ and the streamwise coordinate $x$ have been compressed by a ratio $\kappa = (l_c/h_N)^{2/3}$, where $l_c=\sqrt {\sigma /(\rho g \sin (\beta ))}$ is a capillary length and $h_N$ is the thickness of the uniform film (Nusselt solution). The spatial coordinates $x$ and $y$ have been made dimensionless with the scales $\kappa h_N$ and $h_N$, respectively. The time scale is $\kappa h_N/u_N$, where $u_N=g \sin (\beta ) h_N^2/\nu$ is three times the averaged speed of the Nusselt solution, which corresponds to the speed of kinematic waves at the free surface in the very-long-wave limit. Similarly, the velocity components $u$ and $v$ have been made dimensionless with the scales $u_N$ and $u_N/\kappa$, respectively. This choice of scales enables us to limit significantly the number of coefficients in (2.1) differing from unity. In particular, $\kappa$ is adjusted so that the surface tension term bears a coefficient equal to one.

Three dimensionless groups characterise the flow: the reduced Reynolds number $\delta = 3 Re/\kappa$, where $Re = u_N h_N/(3 \nu )$ is the Reynolds number based on the averaged velocity $u_N/3$; a reduced slope $\zeta = \cot (\beta )/\kappa$; and a parameter $\eta = 1/\kappa ^2$ that compares viscous diffusion and capillary damping of the waves. The reduced Reynolds number $\delta$ compares inertia and viscosity at a scale where surface tension is dominant. Ooshida (Reference Ooshida1999) showed that $\delta$ discriminates between two different regimes for which travelling waves present different tail lengths and properties, which he called the ‘drag–gravity’ and ‘drag–inertia’ regimes. Finally, we introduce the Kapitza number $Ka = \sigma /(\rho \nu ^{4/3} g^{1/3})$, which compares surface tension, viscosity and gravity. This parameter depends only on the fluid properties.

The Nusselt flat-film flow solution can be written as

(2.2)\begin{equation} u = h^2 \left(\frac{y}{h} - \frac{1}{2}\frac{y^2}{h^2} \right) \equiv h^2 g_0({\bar{y}}), \quad \text{where}\ \bar{y} =y/h. \end{equation}

The velocity is then parametrised with only one field, namely the film thickness $h$, which is also the amplitude of the long-wave surface mode of instability. Yet, the velocity distribution of large-amplitude waves may significantly depart from the Nusselt solution. This is an effect of the inertia of the film. An easy way to account for this effect is to parametrise the velocity distribution with another field, usually the flow rate $q=\int _0^h u \,{{\rm d}y}$ or the averaged velocity $\bar {u} = q/h$.

Following Ruyer-Quil et al. (Reference Ruyer-Quil, Kofman, Chasseur and Mergui2014), we decompose the streamwise velocity into a main part and a correction,

(2.3)\begin{equation} u = u^{(0)}(h,q) + \varepsilon u^{(1)}(h,q) \equiv h^2 g_0({\bar{y}}) + \left(\frac{q}{h} - \frac{h^2}{3} \right) f_A(\bar{y}) + \varepsilon u^{(1)}(h,q), \end{equation}

where $f_A$ refers to the ellipse-type velocity profile proposed by Usha et al. (Reference Usha, Chattopadhyay and Tiwari2020):

(2.4)\begin{align} f_A = K_p \left[\sqrt{A^2 - 4(\bar{y} -1)^2} - B \right], \quad\text{with}\ K_p = \frac{1}{\tfrac{1}{4}{A^2} \sin^{{-}1}({2}/{A}) - \tfrac{1}{2} B }, \quad B=\sqrt{A^2 - 4}. \end{align}

The ellipse velocity profile $f_A$ and corrections $u^{(1)}$ fulfil the $O(\varepsilon )$ boundary conditions for the velocity and integral conditions:

(2.5)\begin{equation} \int_0^1 f_A(\bar{y})\, {\rm d}\bar{y} =1,\quad f_A(0)=0,\quad \frac{{\rm d}\, f_A}{{\rm d}{\bar y}}(1)=0, \quad \int_0^h u^{(1)}(h,q)\, {{\rm d}y}=0. \end{equation}

Consequently, (2.3) verifies $q=\int _0^h u\,{{\rm d} y}$ such that this decomposition is unique. Note that $u^{(0)}$ coincides with the Nusselt solution (2.2) for $q=h^3/3$, which corresponds to the in-depth integration of (2.2). The profile $f_A$ involves only one free parameter $A$, which can be varied in the range $]2, \infty [$. In the limit $A\to 2$, $f_A$ approaches a circular shape (see figure 1a). Note that, as $A$ increases, the eccentricity of the ellipse approaches values closer to $1$ (Usha et al. Reference Usha, Chattopadhyay and Tiwari2020), and $f_A$ tends to a parabola as $A\to \infty$. The profiles for higher values of $A\ (= 5, 10)$ are closer to each other and almost coalesce with the parabolic profile.

Figure 1.

(a) Base flow for different values of $A$. (b) Solitary wave profiles corresponding to solutions to EVP and EVPM models (2.8) and (2.9) for $\theta = 90$, $Ka = 3400$, $\delta = 4$ and $A = A^\star =2.23219$.

We retain the non-vanishing terms of highest order, namely first-order inertial terms and second-order viscous diffusion ones. Although second-order inertia terms can be retained, these terms may lead to non-physical behaviour, with the loss of solitary wave solutions, as $\delta$ is raised (Ruyer-Quil et al. Reference Ruyer-Quil, Kofman, Chasseur and Mergui2014):

(2.6a)\begin{gather} \varepsilon u^{(1)}_{yy} = \varepsilon\delta (u^{(0)}_t + u^{(0)}\,u^{(0)}_x + v^{(0)}\,u^{(0)}_y) - b(h) - u^{(0)}_{yy} \nonumber\\ \qquad\, - \,\varepsilon^2 \eta (2u^{(0)}_{xx} + [ u^{(0)}_x|_{y=h} ]_x ) + O(\varepsilon^2), \end{gather}

(2.6b)\begin{gather}\varepsilon u^{(1)}|_{y=0} = 0, \quad \varepsilon u^{(1)}_y|_{y=h} = \varepsilon^2\eta (4 h_x u^{(0)}_x|_{y=h} - v^{(0)}_x|_{y=h} ) + O(\varepsilon^3). \end{gather}

The resulting momentum balance (2.6a) is next averaged with a parabolic weight $w= g_0(\bar {y})$, the solution to $w_{yy} = \mathrm {const.}$, since then

(2.7)\begin{equation} \int_0^h g_0 u^{(1)}_{yy}\, {{\rm d} y} = \tfrac{1}{2} u^{(1)}_y|_{y=h} -\int_0^h u^{(1)}\,{{\rm d}y} = \varepsilon\eta(2h_x u^{(0)}_x|_{y=h} - \tfrac{1}{2} v^{(0)}_x|_{y=h}). \end{equation}

As a result, the computation of the correction $u^{(1)}$, the solution to (2.6), is unnecessary to achieve consistency at order $O(\varepsilon )$. Since $w$ is not proportional to $u^{(0)}$, the WRM differs here from the EIM – see Ruyer-Quil, Chakraborty & Dandapat (Reference Ruyer-Quil, Chakraborty and Dandapat2012) and Samanta, Goyeau & Ruyer-Quil (Reference Samanta, Goyeau and Ruyer-Quil2013) for other examples of the WRM approach where weight and test functions differ. The obtained averaged momentum balance thus reads

(2.8)\begin{gather} \delta q_t = \delta \left[\left(\frac{G}{S} \frac{q^2}{h^2} - \frac{G_1}{S} q h - \frac{G_2}{S} h^4\right)h_x - \left(\frac{F}{S} \frac{q}{h} + \frac{F_1}{S} h^2\right)q_x \right] + \frac{h}{3 S}\left[b(h)- \frac{3q}{h^3}\right] \nonumber\\ \quad + \,\eta \left[\left(\frac{J}{S} \frac{q}{h^2} - \frac{J_1}{S} h\right)h_x^2 - \frac{K}{S} \frac{q_x h_x}{h} - \left(\frac{L}{S} \frac{q}{h} + \frac{L_1}{S} h^2\right) h_{xx} + \frac{M}{S} q_{xx} \right], \end{gather}

where the ordering parameter $\varepsilon$ has been dropped for simplicity. The coefficients $F$, $G$, etc. are functions of the parameter $A$, the expressions of which are given in Appendix A.

Compared to the RQM model derived by Ruyer-Quil & Manneville (Reference Ruyer-Quil and Manneville2000), (2.8) contains deviatory terms that are not present there, corresponding to the coefficients $G_1$, $G_2$, $F_1$, $J_1$ and $L_1$. These terms arise due to the non-self-similar nature of $u_0(h,q)$ in the ansatz (2.3) as the parabolic profile $3g_0$ differs from the elliptic profile $f_A$. In the limit $A\to \infty$, $f_\infty =3g_0$ and these terms vanish. Owing to these deviatory terms, the averaged momentum balance (2.8) differs from the expected form of an averaged momentum balance, akin to a viscous shallow-water equation (Gerbeau & Perthame Reference Gerbeau and Perthame2001), where convective terms are quadratic with respect to the velocity, or equivalently the flow rate $q$, and viscous diffusion terms are linear with respect to $q$.

However, (2.8) can easily be rewritten in a form where convective terms are quadratic and diffusive terms are linear with respect to $q$ by playing with the zeroth-order equivalence $q = h^3/3 + O(\varepsilon )$ to rewrite the inertia and diffusion terms while retaining consistency at $O(\varepsilon )$ for the former and $O(\varepsilon ^2)$ for the latter. The result reads

(2.9)\begin{gather} \delta q_t = \delta \left[\left(\frac{G - 3 G_1 - 9 G_2}{S}\right)\frac{q^2}{h^2} h_x - \left(\frac{F + 3 F_1}{S}\right)\frac{q}{h} q_x \right] + \frac{h}{3 S}\left[b(h)- \frac{3q}{ h^3}\right]\nonumber\\ \qquad\, +\,\eta \left[\left(\frac{J - 3 J_1}{S}\right)\frac{q}{h^2} h_x^2 - \frac{K}{S} \frac{q_x h_x}{h} - \left(\frac{L + 3 L_1}{S} \right) \frac{q}{h} h_{xx} + \frac{M}{S} q_{xx}\right]. \end{gather}

Completed by the (exact) mass balance $h_t + q_x = 0$, (2.8) or (2.9) form two-equation low-dimensional models, referred to as EVP and EVPM models, respectively. These models are consistent up to $O(\varepsilon )$ for inertial terms and up to $O(\varepsilon ^2)$ for diffusive ones. As a consequence, the threshold of instability, $\delta _c = \tfrac {5}{2} \zeta$, is accurately captured by the two models. These two models lead to the same linear stability analysis of the Nusselt flat-film solution (same dispersion relation). Moreover, the locations of the marginal stability curves based on this dispersion relation (not shown) are nearly unaffected by the value of $A$ for a wide range of parameters. This shows that the linear stability of the Nusselt film based on the EVP, EVPM and RQM models are nearly equivalent. Since the RQM model has been shown to accurately capture the solutions to the Orr–Sommerfeld problem (Ruyer-Quil & Manneville Reference Ruyer-Quil and Manneville2002), we conclude that our new formulations reproduce satisfactorily the primary instability of the film for all values of $A>2$. The variations of the coefficients of the driving force ($b$), convective quadratic terms ($q^2 h_x/h$, $q q_x/h$) and linear diffusive terms ($q h_x^2/h^2$, $h_x q_x/h$, $q h_{xx}/h$, $q_{xx}$) with respect to $A$ are illustrated in table 1, which shows an excellent convergence to the RQM model as $A$ is raised.

Table 1.

Coefficients of (2.8) (EVP) and (2.9) (EVPM) compared with the RQM model (equation (14) in Ruyer-Quil et al. Reference Ruyer-Quil, Kofman, Chasseur and Mergui2014).

4. Concluding remarks

Based on a combination of a parabolic and an ellipse profile introduced by Usha et al. (Reference Usha, Chattopadhyay and Tiwari2020), we have derived a new low-dimensional two-equation model for thin-film flows using the WRM. The model is first-order consistent for inertia terms and second-order consistent for diffusive terms, which guarantees that it captures the primary instability adequately. Equation (2.9) presents the usual structure for a shallow-water momentum balance, where convective terms are quadratic with respect to the flow rate $q$ and diffusive terms are linear in $q$. The EVPM (2.9) successfully and satisfactorily captures the nonlinear wave properties of film flows down a vertical wall, in both the drag–gravity and drag–inertia regimes, such as asymptotic wave speed, maximum wave height, wave profiles, streamlines in the moving frame and the flow patterns in the laboratory frame. In particular, it is able to capture the onset of closed separation vortices which can form underneath the troughs of precursory capillary ripples (Dietze et al. Reference Dietze, Al-Sibai and Kneer2009).

This study is another step forward in the modelling of falling thin-film flows, thereby offering new perspectives for theoretical and applied investigations. This calls for further in-depth investigation of the EVPM model and similar formulations, taking full advantage of the adjustable ansatz for the velocity profile.