3.1. Non-dimensionalisation

We first non-dimensionalise the equations according to

(3.1)\begin{equation} \left. \begin{array}{@{}cc@{}} &x = \mathcal{L} \tilde{x},\quad u = \mathcal{U} \tilde{u},\quad H = \epsilon \mathcal{L} \tilde{H},\quad t = \frac{\mathcal{L}}{\mathcal{U}}\tilde{t},\quad y = \epsilon \mathcal{L} \tilde{y},\quad w = \epsilon \mathcal{U} \tilde{w}, \\ &h = \epsilon \mathcal{L} \tilde{h}, \quad p = \frac{\mu \mathcal{U}}{\mathcal{L}}\tilde{p},\quad \sigma = \varSigma +\left (\Delta\varSigma\right ) \tilde{\sigma} , \end{array} \right\} \end{equation}

where $\mathcal {U}$ is some constant characteristic velocity, $\varSigma$ is some constant characteristic surface tension and $\Delta \varSigma$ is a characteristic surface tension variation of interest in the problem. There are three dimensionless parameters: Reynolds number $ Re = {\rho \mathcal {U} \mathcal {L}}/{\mu }$, Marangoni number $\textit {M} = {\Delta \varSigma }/{\epsilon \mu \mathcal {U}}$ and capillary number $\textit {C} = {\epsilon \mu \mathcal {U}}/{\varSigma }$. The parameter $\epsilon$ is included in the definition of $\textit {M}$ and $\textit {C}$ in a way such that (a) the later discussed thresholds for the extensional flow conditions are $\textit {O}(1)$ (see (3.10) and (3.11)) and (b) the complete equations become independent of $\epsilon$ (see § 3.5). Note that in the non-dimensionalisation (3.1), we are considering timescales ${\mathcal {L}}/{\mathcal {U}}$ (see § 4.7 for a discussion of shorter timescales). Henceforth, the tilde $\widetilde {(\ )}$ will be omitted.

3.1.1. Non-dimensional equations

The bulk equations (2.1), (2.2) and (2.3) in dimensionless form are

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

(3.3)$$\begin{gather} Re \left( \frac{\partial u}{\partial t} +u\frac{\partial u}{\partial x}+w\frac{\partial u}{\partial z} \right) ={-} \frac{\partial p}{\partial x}+ \left(\frac{\partial^2 u}{\partial x^2}+\frac{1}{\epsilon^2}\frac{\partial^2 u}{\partial z^2}\right), \end{gather}$$

(3.4)$$\begin{gather} Re \left( \frac{\partial w}{\partial t} +u\frac{\partial w}{\partial x}+w\frac{\partial w}{\partial z} \right) ={-} \frac{1}{\epsilon^2}\frac{\partial p}{\partial z}+ \left(\frac{\partial^2 w}{\partial x^2}+\frac{1}{\epsilon^2}\frac{\partial^2 w}{\partial z^2}\right). \end{gather}$$

Next, the kinematic boundary conditions (2.4) are given by

(3.5)\begin{equation} w|_{z=z^{{\pm}}}=\left(\frac{\partial H}{\partial t}\pm\frac{1}{2}\frac{\partial h}{\partial t}\right) + u|_{z=z^{{\pm}}}\left(\frac{\partial H}{\partial x}\pm\frac{1}{2}\frac{\partial h}{\partial x}\right). \end{equation}

With the non-dimensionalisation, the normal and tangential stress boundary conditions (2.5) and (2.6) are given by

(3.6) $$\begin{gather} \frac{\epsilon^2}{C}\left(1+MC\sigma_{{\pm}}\right)\kappa_{{\pm}} = \left(p-\frac{2\left(\dfrac{\partial H}{\partial x}\pm\dfrac{1}{2}\dfrac{\partial h}{\partial x}\right)^2}{1+\epsilon^2\left(\dfrac{\partial H}{\partial x}\pm\dfrac{1}{2}\dfrac{\partial h}{\partial x}\right)^2}\epsilon^2\frac{\partial u}{\partial x} \right.\nonumber\\ +\left.\left.\frac{2\left(\dfrac{\partial H}{\partial x}\pm\dfrac{1}{2}\dfrac{\partial h}{\partial x}\right)}{1+\epsilon^2\left(\dfrac{\partial H}{\partial x}\pm\dfrac{1}{2}\dfrac{\partial h}{\partial x}\right)^2}\left(\frac{\partial u}{\partial z}+\epsilon^2\frac{\partial w}{\partial x}\right)-\frac{2}{1+\epsilon^2\left(\dfrac{\partial H}{\partial x}\pm\dfrac{1}{2}\dfrac{\partial h}{\partial x}\right)^2}\frac{\partial w}{\partial z}\right)\right|_{z=z^{{\pm}}} \end{gather}$$

and

(3.7)\begin{align} \epsilon^2M \frac{\partial \sigma_{{\pm}}}{\partial x}&=\left({\pm}\frac{2\left(\dfrac{\partial H}{\partial x}\pm\dfrac{1}{2}\dfrac{\partial h}{\partial x}\right)\epsilon^2}{\sqrt{1+\epsilon^2\left(\dfrac{\partial H}{\partial x}\pm\dfrac{1}{2}\dfrac{\partial h}{\partial x}\right)^2}}\left(\frac{\partial w}{\partial z}-\frac{\partial u}{\partial x}\right) \right.\nonumber\\ &\quad \pm\left.\left.\frac{1-\epsilon^2\left(\dfrac{\partial H}{\partial x}\pm\dfrac{1}{2}\dfrac{\partial h}{\partial x}\right)^2}{\sqrt{1+\epsilon^2\left(\dfrac{\partial H}{\partial x}\pm\dfrac{1}{2}\dfrac{\partial h}{\partial x}\right)^2}}\left(\frac{\partial u}{\partial z}+\epsilon^2\frac{\partial w}{\partial x}\right)\right)\right|_{z=z^{{\pm}}}, \end{align}

where the curvatures $\kappa _{\pm }$ (2.7) are given by

(3.8)\begin{equation} \kappa_{{\pm}} = \frac{\mp\dfrac{\partial^2 H}{\partial x^2}-\dfrac{1}{2}\dfrac{\partial^2 h}{\partial x^2}}{\left(1+\epsilon^2\left(\dfrac{\partial H}{\partial x}\pm\dfrac{1}{2}\dfrac{\partial h}{\partial x}\right)^2\right)^{{3}/{2}}}. \end{equation}

These equations involve four non-dimensional parameters, $\epsilon$, $ Re $, $M$ and $C$.

3.2. Conditions for extensional flow

In order to have extensional flow, i.e. $u(x,z,t) \approx u(x,t)$ and $p(x,z,t) \approx p(x,t)$ to leading order in $\epsilon$, we require three conditions. The first condition is that inertia does not play a dominant role:

(3.9)\begin{equation} Re \left( \frac{\partial u}{\partial t} +u\frac{\partial u}{\partial x}+w\frac{\partial u}{\partial z} \right)= \textit{O}(1) \text{ or less}. \end{equation}

Starting with (3.3), this condition yields ${\partial ^2 u}/{\partial z^2} = 0$ to leading order; see (3.13). The second condition is that the tangential stress is not too large:

(3.10)\begin{equation} M \frac{\partial \sigma_{{\pm}}}{\partial x} = \textit{O}(1) \text{ or less}. \end{equation}

This condition yields the horizontal velocity field $u(x,z,t) \approx u(x,t)$; see (3.14). The third condition is that the capillary pressure due to surface tension is not too large:

(3.11)\begin{equation} \frac{1}{C}\left(1+MC\sigma_\pm\right)\kappa_{{\pm}} = \textit{O}(1) \text{ or less}. \end{equation}

This condition yields the pressure field $p(x,z,t) \approx p(x,t)$; see (3.19).

The first two conditions (3.9) and (3.10) are commonly discussed elsewhere (De Wit et al. Reference De Wit, Gallez and Christov1994; Breward Reference Breward1999; Timmermans & Lister Reference Timmermans and Lister2002; Bowen & Tilley Reference Bowen and Tilley2013; Kitavtsev et al. Reference Kitavtsev, Fontelos and Eggers2018). The third condition (3.11) has also been mentioned (Howell Reference Howell1996; Breward Reference Breward1999). In the case of symmetry ($\sigma = \sigma _{\pm }$, $\kappa = \kappa _{\pm }$), the right-hand side of (3.11) may be replaced with the phrase ‘$\textit {O}(\epsilon ^{-2})$ or less’ (see Appendix A). Thus, (3.11) is a condition that becomes more strict due to the asymmetry of interest.

3.3. Leading-order expansion

We expand the horizontal velocity as

(3.12)\begin{equation} u(x,z,t) = u_0(x,z,t) + \epsilon^2 u_1(x,z,t) + \cdots \end{equation}

and consider analogous expansions for $w, H, h,$ and $p$. We consider the leading-order expansion first where we ignore relative errors of order $\textit {O}(\epsilon ^2)$. The horizontal momentum equation (3.3) using the first condition (3.9) gives

(3.13)\begin{equation} \frac{\partial^2 u_{0}}{\partial z^2} = 0. \end{equation}

The tangential boundary conditions (3.7), after using the second condition (3.10), give

(3.14)\begin{equation} \left.\frac{\partial u_0}{\partial z}\right|_{z=z^{{\pm}}} = 0,\end{equation}

which shows that the leading-order horizontal velocity is independent of $z$,

(3.15)\begin{equation} u_0 = u_0(x,t). \end{equation}

Then, continuity (3.2) gives

(3.16)\begin{equation} w_0 = \bar{w}(x,t) - (z-H_0)\frac{\partial u_0}{\partial x}, \end{equation}

where $\bar {w}$ denotes the average of $w$ in the vertical direction $(\,\bar{f}:= ({1}/{h})\int _{z^-}^{z^+} f(x,z,t)\,\textrm {d}z)$.

Next, turning to the vertical momentum equation (3.4) and using (3.9) gives

(3.17)\begin{equation} \frac{\partial p_0}{\partial z} = \frac{\partial^2 w_0}{\partial z^2}, \end{equation}

which upon substitution of (3.16) yields

(3.18)\begin{equation} p_0 = p_0(x,t). \end{equation}

Then, considering the normal stress boundary conditions (3.6), along with continuity (3.2), gives

(3.19)\begin{equation} 0 = p_0 + 2 \frac{\partial u_0}{\partial x}. \end{equation}

Finally, with the kinematic boundary conditions (3.5), and using (3.16) for $w_0$, we find

(3.20)\begin{equation} \bar{w} \mp \frac{1}{2}h_0 u_0 = \left(\frac{\partial H_0}{\partial t} \pm \frac{1}{2}\frac{\partial h_0}{\partial t}\right) + u_0 \left(\frac{\partial H_0}{\partial x} \pm \frac{1}{2}\frac{\partial h_0}{\partial x}\right), \end{equation}

which can be added and subtracted to find

(3.21)\begin{equation} \frac{\partial H_0}{\partial t} + u_0 \frac{\partial H_0}{\partial x} = \bar{w} \end{equation}

and

(3.22)\begin{equation} \frac{\partial h_0}{\partial t} + \frac{\partial}{\partial x}\left(u_0h_0\right) = 0.\end{equation}

3.4. Next-order expansion

We can follow the same steps as in § 3.3 for the $\textit {O}(\epsilon ^2)$ equations. The details are provided in the supplementary material. Explicitly, we consider the horizontal momentum equation and the tangential stress conditions to deduce

(3.23)\begin{equation} Re h_0\left(\frac{\partial u_0}{\partial t}+u_0 \frac{\partial u_0}{\partial x}\right) = M\left(\frac{\partial \sigma_+}{\partial x} + \frac{\partial \sigma_-}{\partial x}\right) + 4 \frac{\partial}{\partial x}\left(h_0\frac{\partial u_0}{\partial x} \right).\end{equation}

Similarly, the vertical momentum equation and the normal stress conditions lead to

(3.24)$$\begin{gather} Re h_0 \left(\frac{\partial \bar{w}}{\partial t} + u_0\frac{\partial \bar{w}}{\partial x}\right) =\frac{1}{2}M\frac{\partial^2 }{\partial x^2}\left(h_0(\sigma_+{-}\sigma_-)\right)+\frac{\partial}{\partial x}\left(\frac{\partial H_0}{\partial x}\left(\frac{2}{C}+M\left(\sigma_+{+}\sigma_-\right)\right)\right) \nonumber\\ +\, 4 \frac{\partial}{\partial x}\left(h_0 \frac{\partial H_0}{\partial x}\frac{\partial u_0}{\partial x}\right). \end{gather}$$

3.5. Complete equations

In this subsection, the complete non-dimensional equations are given for clarity. From this point onwards, we refer to $h_0$ as $h$, $H_0$ as $H$ and $u_0$ as $\bar {u}$. Then, the evolution equations for the dimensional counterparts $h$, $H$, $\bar {u}$ and $\bar {w}$ are given by

(3.25)$$\begin{gather} \frac{\partial h}{\partial t} + \frac{\partial}{\partial x}\left(\bar{u}h\right) = 0, \end{gather}$$

(3.26)$$\begin{gather}\frac{\partial H}{\partial t} + \bar{u} \frac{\partial H}{\partial x} = \bar{w}, \end{gather}$$

(3.27)$$\begin{gather} Re h\left(\frac{\partial \bar{u}}{\partial t}+\bar{u} \frac{\partial \bar{u}}{\partial x}\right) = M\left(\frac{\partial \sigma_+}{\partial x} + \frac{\partial \sigma_-}{\partial x}\right) + 4 \frac{\partial}{\partial x}\left(h\frac{\partial \bar{u}}{\partial x} \right), \end{gather}$$

(3.28)$$\begin{gather} Re h \left(\frac{\partial \bar{w}}{\partial t} + \bar{u}\frac{\partial \bar{w}}{\partial x}\right) =\frac{1}{2}M\frac{\partial^2 }{\partial x^2}\left(h(\sigma_{+} - \sigma_-)\right)+\frac{\partial}{\partial x}\left(\frac{\partial H}{\partial x}\left(\frac{2}{C}+M\left(\sigma_{+} +\sigma_{-}\right)\right)\right) \nonumber\\ +\, 4 \frac{\partial}{\partial x}\left(h \frac{\partial H}{\partial x}\frac{\partial \bar{u}}{\partial x}\right). \end{gather}$$

In the above equations, the dynamics depend only on $ Re $, $M$ and $C$. In the symmetric case $H=0$, $\sigma _{\pm } = \sigma$ and $\bar {w} = 0$, we recover equations already familiar in the literature (De Wit et al. Reference De Wit, Gallez and Christov1994). The capillary pressure term $\tfrac {1}{2}({\epsilon ^2}/{C}) h ({\partial ^3 h}/{\partial x^3})$ does not appear in the above equations since we consider $C = \textit {O}(1)$ or greater; see Appendix A.

Equations (3.25) and (3.26) have been given in the literature previously (Howell Reference Howell1996; Breward Reference Breward1999). In addition, (3.27) was given by Breward (Reference Breward1999), though inertia was omitted, and (3.28) with $\sigma _+ = \sigma _-$ constant was given by Howell (Reference Howell1996). To the best of the authors’ knowledge, the evolution equation (3.28) for $\bar {w}$ with Marangoni effects has not been considered so far in the literature. In the example of surfactant deposition on a sheet, which we consider in §§ 4 and 5, asymmetric surface tension variations are considered in (3.28) in order to solve for the resulting centreline deformation $H$.

Before proceeding further, we provide intuition for (3.25)–(3.28). The evolution equation for $h$ (3.25) expresses conservation of mass. The evolution equation for $H$ (3.26) can be regarded as ${\textrm {D}H}/{\textrm {D}t}= \bar {w}$ where ${\textrm {D}}/{\textrm {D}t}$ is the material derivative. Thus, the evolution equation for $H$ indicates that the centreline moves with the mean vertical flow. For the horizontal momentum equation (3.27), the left-hand side is inertia, the first term of the right-hand side is from the Marangoni force and the second term of the right-hand side is viscous stresses (sometimes referred to as ‘Trouton viscosity’). The vertical momentum equation (3.28) can be discussed in a similar way.

4.1. Non-dimensional equations

Since there is no background flow, we note that $\mathcal {U}$ is set by the tangential Marangoni stress. Thus, we take $\mathcal {U} = {\Delta \varSigma }/{\epsilon \mu }$ (in other words, take $M =1$). Then, $ Re = {\rho \mathcal {L}\Delta \varSigma }/{\epsilon \mu ^2}$, ${C = {\Delta \varSigma }/{\varSigma }}$. We define $\epsilon = {h_I}/{\mathcal {L}}$. In addition, non-dimensionalise the surfactant concentration with a characteristic surfactant concentration $\varGamma _c$. Explicitly,

(4.1)\begin{equation} \left. \begin{array}{@{}cc@{}} &x = \mathcal{L} \tilde{x}, \quad u = \dfrac{\Delta \varSigma}{\epsilon \mu}\tilde{u},\quad H = \epsilon \mathcal{L} \tilde{H},\quad t = \dfrac{\epsilon \mu \mathcal{L}}{\Delta \varSigma}\tilde{t},\quad y = \epsilon \mathcal{L} \tilde{y},\quad w = \dfrac{\Delta \varSigma}{\mu} \tilde{w}, \\ &h = \epsilon \mathcal{L} \tilde{h}, \quad \sigma = \varSigma+ \Delta \varSigma \tilde{\sigma}, \quad \varGamma = \varGamma_c \tilde{\varGamma}. \end{array} \right\} \end{equation}

Note that, by definition, $\tilde {\sigma }(\varGamma = 0)=0$. Again, the tilde $\widetilde {(\ )}$ is dropped. We consider $C = \textit {O}(1)$ or larger so that we have extensional flow (see § 3.2). The horizontal length scale $\mathcal {L}$ is given by the width of variation of the initial surfactant concentration. For example a Gaussian initial surfactant distribution $\varGamma _I$ is represented as $\varGamma _I \propto \textrm {e}^{-x^2}$.

Equations (3.25) and (3.26) are as before. Since $ Re \ll 1$, the inertia terms of (3.27) and (3.28) can be neglected to deduce, respectively, quasi-static balances of horizontal momentum and vertical momentum given by

(4.2)$$\begin{gather} 0= \left(\frac{\partial \sigma_+}{\partial x} + \frac{\partial \sigma_-}{\partial x}\right) + 4 \frac{\partial}{\partial x}\left(h\frac{\partial \bar{u}}{\partial x} \right), \end{gather}$$

(4.3)$$\begin{gather} 0=\frac{1}{2}\frac{\partial^2 }{\partial x^2}\left(h(\sigma_{+}-\sigma_{-})\right)+\frac{\partial}{\partial x}\left(\frac{\partial H}{\partial x}\left(\frac{2}{C}+\left(\sigma_{+}+\sigma_-\right)\right)\right) \nonumber\\+\ 4 \frac{\partial}{\partial x}\left(h \frac{\partial H}{\partial x}\frac{\partial \bar{u}}{\partial x}\right). \end{gather}$$

It should be noted that there is now a single parameter ($C$) remaining, which as noted at the beginning of the section, only changes the prefactor of $H$ and $\bar {w}$. The quasi-static balance occurs on the dimensional timescale given by $t_{qs}:={\epsilon \mu \mathcal {L}}/{\Delta \varSigma }$. The system adjusts itself from the initial condition to the quasi-static balance on dimensional timescales much shorter than $t_{qs}$ and these details are discussed in § 4.7, where inertia is important. Since diffusion of surfactants is neglected, the equations for surfactant concentration at the top and bottom surfaces, $\varGamma _{\pm }(x,t)$, are given to leading order, correct up to $\textit {O}(\epsilon ^2)$ errors, by

(4.4)\begin{equation} \frac{\partial \varGamma_{{\pm}}}{\partial t} + \frac{\partial}{\partial x}\left(\bar{u}\varGamma_{{\pm}}\right)=0.\end{equation}

If diffusion is not ignored, the term $({\epsilon \mu D}/{ \mathcal {L} \Delta \varSigma }) ({\partial ^2 \varGamma _{\pm }}/{\partial x^2})$ (diffusion coefficient $D$) appears on the right-hand side of (4.4). In neglecting diffusion of surfactants, we are assuming that the Péclet number ${\mathcal {L}\Delta \varSigma }/{\epsilon \mu D}$ is sufficiently large. Neglecting the influence of diffusion is also physically motivated. Upon taking surfactant diffusivity $D = 10^{-10}$ m$^2$ s$^{-1}$, along with the values quoted at the beginning of the section ($\Delta \varSigma = 0.006$ N m$^{-1}$, $\mathcal {L} = 100\,\mathrm {\mu } \textrm {m}$, $\mu = 1$ Pa s, $\epsilon = 0.1$) and the Péclet number ${\mathcal {L} \Delta \varSigma }/{\epsilon \mu D} \approx 6 \times 10^4 \gg 1$.

4.2. Initial and boundary conditions

For the initial conditions, we consider a uniform sheet otherwise at rest, along with some initial surfactant distribution $\varGamma _{I\pm }(x)$ where $\varGamma _{I\pm } \rightarrow 0$ as $x \rightarrow \pm \infty$. Far field conditions are taken for the boundary condition. See figure 3 for a particular example where $\varGamma _{I+}=\textrm {e}^{-x^2}$ and $\varGamma _{I-} = 0$, which is considered in § 4.4. Explicitly,

(4.5)\begin{equation} h = 1,\ H = 0,\ \bar{u} = 0,\ \bar{w} = 0,\ \varGamma_{{\pm}} = \varGamma_{I\pm} \quad \text{at } t = 0\end{equation}

and

(4.6)\begin{equation} h = 1,\ H = 0,\ \frac{\partial \bar{u}}{\partial x}= 0,\ \bar{w} = 0,\ \varGamma_{{\pm}} = 0 \quad \text{at } x ={\pm} \infty. \end{equation}

Figure 3. Initial set-up considered in § 4.4. An initial Gaussian surfactant distribution $\varGamma _{I+}(x) = \textrm {e}^{-x^2}$, $\varGamma _{I-}(x) = 0$ is deposited onto the top of a uniform sheet at rest $-\tfrac {1}{2} \leq z \leq \tfrac {1}{2}$ at $t = 0$. The horizontal direction is given by the $x$ axis and the vertical direction is given by the $z$ axis. See (4.1) for the non-dimensionalisation.

4.3. Lagrangian coordinates

Lagrangian coordinates are useful in discussing the inertialess evolution of liquid sheets and liquid threads (Eggers & Fontelos Reference Eggers and Fontelos2015). In this problem, using Lagrangian coordinates allows us to deduce physical conclusions for a general surfactant isotherm $\sigma = \sigma (\varGamma )$.

Before switching to Lagrangian coordinates, we first do some algebra. Integrating the horizontal momentum equation (4.2) with respect to $x$, and using boundary conditions (4.6), gives

(4.7)\begin{equation} \frac{\partial \bar{u}}{\partial x} ={-}\frac{1}{4h(x,t)} \left(\sigma(\varGamma_+(x,t))+\sigma(\varGamma_-(x,t))\right),\end{equation}

which can be substituted into (4.3) and integrated twice with respect to $x$ to deduce an analytic expression for the centreline $H$ shape given by

(4.8)\begin{equation} H(x,t) ={-}\frac{C}{4}h(\sigma(\varGamma_+(x,t))-\sigma(\varGamma_-(x,t))). \end{equation}

Denote the Lagrangian coordinates by the initial material spatial coordinate $s$ defined so that $x=x(s,t)$ with $s = x(s,0)$. At initial time $t = 0$, the Eulerian and Lagrangian coordinates therefore agree. Denoting the Lagrangian time derivative by ${\textrm {D}}/{\textrm {D}t}$, (3.25) and (4.4) give

(4.9)\begin{equation} \left.\frac{{\rm D}h}{{\rm D}t}\right|_{(s,t)} = \left.- h \frac{\partial \bar{u}}{\partial x} \right|_{(x(s,t),t)}\end{equation}

and

(4.10)\begin{equation} \left.\frac{{\rm D} \varGamma_{{\pm}}}{{\rm D}t}\right|_{(s,t)} ={-} \left.\varGamma_{{\pm}}\frac{\partial \bar{u}}{\partial x}\right|_{(x(s,t),t)}, \end{equation}

which gives that

(4.11)\begin{equation} \left.\frac{{\rm D}}{{\rm D}t}\left(\frac{\varGamma_{{\pm}}}{h}\right)\right|_{(s,t)} = \left.\frac{1}{h^2}\left(- \varGamma_{{\pm}}h \frac{\partial \bar{u}}{\partial x}+ h \varGamma_{{\pm}} \frac{\partial \bar{u}}{\partial x}\right)\right|_{(x(s,t),t)} = 0.\end{equation}

Thus, in Lagrangian coordinates,

(4.12)\begin{equation} \varGamma_{{\pm}}(s,t) = \varGamma_{I\pm}(s)h(s,t);\end{equation}

note that (4.12) would not have been deduced if diffusion was included. The result makes sense physically by considering a material volume of initial width $\Delta s$. For example, halving $h$ would mean that the width of the volume has doubled, the interfacial area along a surface doubles, and thus in the case of no diffusion, surfactant concentration would be halved. The conservation equation (4.11) has been discussed by Chomaz (Reference Chomaz2001). Finally, (4.7) and (4.12) can be substituted into (4.9) to deduce that

(4.13)\begin{equation} \left.\frac{{\rm D}h}{{\rm D}t}\right|_{(s,t)} =\frac{1}{4}\left(\sigma(\varGamma_{I+}(s)h(s,t))+\sigma(\varGamma_{I-}(s)h(s,t))\right) .\end{equation}

Equation (4.13) is a first-order ordinary differential equation (ODE) for $h$ in Lagrangian coordinates at a given point $s$ with the initial condition $h=1$. An expression for $\bar {w}$ can also be found in terms of $h$ from

(4.14)\begin{equation} \bar{w}(s,t) = \left.\frac{{\rm D}H}{{\rm D}t}\right|_{(s,t)} , \end{equation}

which directly follows from (3.26). Thus, given a surfactant isotherm $\sigma = \sigma (\varGamma )$, one only needs to solve a single ODE, (4.13), to find the evolution of $h$ following a given Lagrangian coordinate $(s,t)$, which in turn can be used to find the evolution of $H, \bar {u}, \bar {w}$, and $\varGamma _{\pm }$ via, respectively, (4.8), (4.7), (4.14) and (4.12). Finally, the Lagrangian and Eulerian coordinates are related via

(4.15)\begin{equation} \left.\frac{\partial x}{\partial s}\right|_{(s,t)} = \frac{1}{h(s,t)}, \end{equation}

which can be deduced from standard arguments using conservation of mass (see Appendix D).

4.4. Example: linear isotherm, Gaussian deposition on top

In this subsection, we consider an example where the initial surfactant distribution given by $\varGamma _{I+}(s) = \textrm {e}^{-s^2}$ and $\varGamma _{I-}(s) = 0$ and the surface tension isotherm is linear ${\sigma (\varGamma ) = - \varGamma}$ (see figure 3 for the set-up). Then, we have

(4.16)\begin{equation} \left.\frac{{\rm D}h}{{\rm D}t}\right|_{(s,t)} ={-} K(s) h(s,t), \end{equation}

where $K(s):= \tfrac {1}{4}\textrm {e}^{-s^2}$. Using the initial condition $h(s,0) = 1$, we have the analytic solution

(4.17)\begin{equation} h(s,t) = {\rm e}^{{-}K(s)t}. \end{equation}

Equation (4.17) is the solution for $h$ given the (3.25), (3.26), (3.27), (3.28), (4.5) and (4.6) with $\sigma (\varGamma ) = - \varGamma$, $\varGamma _{I+}(s)=\textrm {e}^{-s^2}$ and $\varGamma _{I-}(s) = 0$. In turn, the solution for $h$ given by (4.17) can be substituted into (4.8), (4.7), (4.14) and (4.12) to deduce equations for $H$, $\bar {u}$, $\bar {w}$ and $\varGamma _{\pm }$ via

(4.18)\begin{gather} H(s,t) = CK(s) {\rm e}^{{-}2K(s)t}, \end{gather}

(4.19)\begin{gather} \bar{u}(s,t) = \int_0^{s} K(s'){\rm e}^{K(s')t}\,{\rm d}s', \end{gather}

(4.20)\begin{gather} \bar{w}(s,t) ={-}2 C K^2(s) {\rm e}^{{-}2K(s)t}, \end{gather}

(4.21a,b)\begin{gather} \varGamma_+(s,t) = 4K(s){\rm e}^{{-}K(s)t}\quad \text{and}\quad \varGamma_-(s,t) =0. \end{gather}

The corresponding transformation into Eulerian coordinates is given by

(4.22)\begin{equation} x(s,t) = \int_0^{s}\frac{1}{h(s',t)} \,{\rm d}s', \end{equation}

which follows from (4.15) and the symmetry condition $x(0,t)=0$.

In particular, we can conclude from (4.13) that $h$ does not pinch off in finite time. A point of caution is that as the thickness of the sheet tends towards rupture, other forces such as van der Waals forces become important. The rupture of a liquid film with van der Waals forces has been discussed elsewhere, both for the case without surfactants (Vaynblat, Lister & Witelski Reference Vaynblat, Lister and Witelski2001) and with surfactants (Wee et al. Reference Wee, Wagoner and Basaran2022). In comparison, when the initial condition is a uniformly thinning sheet (see § 5.4), there is finite time pinch off even without van der Waals forces.

4.4.1. Description of the solution

Next, we give a description of the solution given by (4.17), (4.18), (4.19), (4.20) and (4.21a,b) in Eulerian coordinates to help gain intuition for the thinning and rupture of the sheet. For the figures in this section, we pick a particular choice of the capillary number, $C = 0.5$. Other choices for $C$ with $C = \textit {O}(1)$ or greater would have been valid also. The effect of changing $C$ is only a change in the prefactor of $H$ and $\bar {w}$ (see (4.8) and (4.14)).

An example with $C = 0.5$ is given in figure 4. In particular, figure 4(a) shows the bottom of the sheet $H-\tfrac {1}{2}h$, the centreline $H$, the top of the sheet $H+ \tfrac {1}{2}h$ and the surfactant distribution $\varGamma _{+}$ at the initial time. By initial time, we mean $t = 0$ where a quasi-static balance has already been achieved and so, in particular, the centreline is bent (see §§ 4.1 and 4.7). Figure 4(b) shows the velocity profiles $\bar {u}$ and $\bar {w}$ at the initial time. Similarly, figure 4(c,d) report the variables at $t = 1$ and figure 4(e, f) report the variables at $t = 10$. The times $t = 1$ and $10$ were chosen to represent mid and late times, respectively, for the thinning problem.

Figure 4. Plots of the space and time evolution of the thin-film dynamics, as given by (4.17), (4.22), (4.18), (4.19), (4.20) and (4.21a,b) for $C=0.5$. The spatial variable $x$ is in Eulerian coordinates. The three curves from bottom to top in (a) are $H - \tfrac {1}{2}h$ (black), $H$ (blue) and $H + \tfrac {1}{2}h$ (coloured with surfactant concentration $\varGamma$) at the initial time and (b) shows the velocity profiles $\bar {u}$ (blue) and $\bar {w}$ (red) at the initial time. Similarly, (c,d) report the variables at $t = 1$ and (e, f) report the variables at $t = 10$. For the non-dimensionalisation, see (4.1).

The non-uniform surfactant distribution means that there are Marangoni effects, where surfactants spread from a region of high concentration to low concentration (see figure 4a,c,e), since low concentration regions have higher surface tension. As result of the horizontal velocity $\bar {u}$ induced by the spreading of surfactants, $h$ decreases and the sheet thins as seen in figure 4(a,c,e).

Since asymmetry has been accounted for, it is also possible to discuss the evolution of the centreline $H$ and the vertical velocity $\bar {w}$. There is only non-uniform surface tension on the top of the sheet. Initially, the centreline $H$ bends towards the top surface, which makes sense physically since the surface tension on the top surface is lower (the same fact can be seen from (4.8)). For the example shown in figure 4, the top of the sheet is pulled away from $x = 0$. As result, there is a region of negative $\bar {w}$ around $x = 0$. Then, as time progresses, $H$ reverts back towards the straight configuration at $z = 0$.

4.5. Example with a satellite drop: linear isotherm, double maximum

The example shown in § 4.4 considered a single maximum in the initial surfactant distribution. Another important example would be to consider multiple maxima. Continuing from § 4.4, we consider an example where the initial surfactant distribution is given by $\varGamma _{I+}(s) = \textrm {e}^{-(s-l)^2}+\textrm {e}^{-(s+l)^2}$ for some $l>0$ and $\varGamma _{I-}(s) = 0$, which for $l\gtrsim 0.71$ has two maxima at $\pm s_m$ with $s_m >0$ and so we expect pinching to occur at two points $s_r = \pm s_m$ (see § 4.6.1), with the consequence being that ‘satellite drops’ are formed. The parameter $l$ is thus a geometric variable that controls how far apart the maxima are located. Upon replacing the definition of $K(s)$ with $K(s):=\tfrac {1}{4}(\textrm {e}^{-(s-l)^2}+\textrm {e}^{-(s+l)^2})$, expressions for $h, H, \bar {u}, \bar {w},$ and $\varGamma _{\pm }$ are still given, respectively, by (4.17), (4.18), (4.19), (4.20) and (4.21a,b) along with conversion into Eulerian coordinates via (4.22). In order to discuss a real drop, one should consider the axisymmetric case rather than the two-dimensional set-up that is considered in this paper. As mentioned in the introduction, the formation of satellite drops due to rupture via Marangoni effects is an interesting avenue of research.

An example with parameters $C = 0.5$ and $l = 2$ is shown in figure 5. In particular, figure 5(a) shows the initial surfactant concentration $\varGamma _{I+}(s) = \textrm {e}^{-(s-2)^2}+\textrm {e}^{-(s+2)^2}$. From (4.17), we have that the point of minimum thickness (pinch off) rupture occurs near $s_m \approx \pm 2$. For example, figure 5(b) shows the sheet profile $H \pm \tfrac {1}{2}h$ (black curves) at $t = 10$, with a satellite drop around $x = 0$.

Figure 5. A deposition problem with initial surfactant concentration given by $\varGamma _{I+}(x) = \textrm {e}^{-(x-2)^2}+\textrm {e}^{-(x+2)^2}$, $\varGamma _{I-}(x) = 0$. (a) The initial surfactant concentration at the top surface $\varGamma _{I+}$. (b) The sheet profile $H \pm \tfrac {1}{2}h$ at $t = 20$, where the inset shows the sheet profile $H\pm \tfrac {1}{2}h$ near the satellite drop (included to illustrate the drop shape). In (b), the inset and the main figure have the same axes. The spatial variable $x$ is in Eulerian coordinates. The isotherm chosen is $\sigma (\varGamma ) = - \varGamma$ and the capillary parameter $C = 0.5$. For the non-dimensionalisation, see (4.1). A ‘satellite drop’ is formed around $x = 0$.

A point of physical interest in areas such as aerosol production at the ocean–air interface is the volume of the satellite droplet that is formed. In our framework, the volume (since we have a two-dimensional set-up, really, a surface area) of the droplet produced, $V$, can be written as

(4.23)\begin{equation} V = \int_{x({-}s_m,t_r)}^{x(s_m,t_r)} h(x,t_r)\, {{\rm d}\kern6em x} = \int_{{-}s_m}^{s_m} \,{\rm d}s = 2s_m. \end{equation}

For $l\gg 1$, $s_m \approx l$ and then $V \approx 2l$ and the volume of the satellite drop produced is just given by the spacing between the maxima (in fact, $l$ does not need to be much greater than $1$, e.g. $l=1.4$ already gives $s_m = 1.399$ to 3 decimal places). In general, the same argument gives that if pinching occurs at two locations with Lagrangian coordinates $s_L< s_R$, then the volume of satellite drop produced is $V = s_R - s_L$.

4.6. Comparison between symmetry and asymmetry

In this section, we investigate the differences and similarities of a top–bottom asymmetric configuration and a top–bottom symmetric configuration for general surface tension isotherms $\sigma = \sigma (\varGamma )$. More explicitly, we compare the thinning and rupture of two cases given by the purely asymmetric case (PA),

(4.24a,b)\begin{equation} \varGamma_{I+}^{PA}(s) = f(s) \quad\text{and}\quad \varGamma_{I-}^{PA}(s)=0, \end{equation}

and the ‘corresponding’ symmetric case (S),

(4.25a,b)\begin{equation} \varGamma_{I+}^S(s) = \frac{1}{2}f(s)\quad\text{and}\quad \varGamma_{I-}^S(s)=\frac{1}{2}f(s). \end{equation}

A function $F(x)$ is said to be strictly monotonically decreasing if $a>b \Rightarrow F(a) < F(b)$. In this section, we only consider surfactant isotherms $\sigma = \sigma (\varGamma )$ that are strictly monotonically decreasing. In other words, the addition of any amount of surfactant strictly decreases the surface tension. In § 4.6.1, we show that the rupture occurs at the Lagrangian point where $f$ is maximal for both the purely asymmetric case (PA) and the symmetric case (S). In § 4.6.2, we show that the purely asymmetric case (PA) thins slower/faster than the symmetric case (S) when $\sigma = \sigma (\varGamma )$ is convex/concave (see figure 6). For comments about cases where there are surfactants on each side but not necessarily deposited symmetrically, see Appendix E, where it is found that the conclusions about the rate of thinning and convexity/concavity still hold. Experimentally, convexity and concavity can be observed for various surface tension isotherms (Prosser & Franses Reference Prosser and Franses2001; Liu & Duncan Reference Liu and Duncan2003; Erinin et al. Reference Erinin, Liu, Liu, Mostert, Deike and Duncan2023).

Figure 6. Schematic diagram of an isotherm $\sigma = \sigma (\varGamma )$ that is (a) convex or (b) concave. In general, we consider $\sigma = \sigma (\varGamma )$ strictly monotonically decreasing.

4.6.1. Location of minimum thickness

Here, we analyse the location of minimum thickness of the sheet $h$ at a given time. In particular, this will give the location at which there is pinching. From (4.13), we have

(4.26)\begin{equation} \left.\frac{{\rm D}h^{PA}}{Dt}\right|_{(s,t)} = \frac{1}{4}\sigma(f(s)h^{PA}(s,t)) \end{equation}

and

(4.27)\begin{equation} \left.\frac{{\rm D}h^S}{{\rm D}t}\right|_{(s,t)} = \frac{1}{2}\sigma\left(\frac{1}{2}f(s)h^S(s,t)\right). \end{equation}

We first analyse the purely asymmetric case (PA). Consider two Lagrangian points denoted by $s_1$ and $s_2$ with $f(s_1)>f(s_2)$. Suppose that at some $t = t_e$, the thickness is equal at the two Lagrangian points $h^{PA}(s_1,t_e) = h^{PA}(s_2,t_e)$. Then, it follows from (4.26) that

(4.28)$$\begin{gather} \frac{{\rm D}}{{\rm D}t}\left.\left(h^{PA}(s_1,t)-h^{PA}(s_2,t)\right)\right|_{t = t_e} = \frac{1}{4}\left(\sigma\left(f(s_1)h^{PA}(s_1,t_e)\right) - \sigma\left(f(s_2)h^{PA}(s_2,t_e)\right)\right)\nonumber\\ < 0, \end{gather}$$

where the inequality follows since $\sigma$ is strictly monotonically decreasing and $f(s_1)>f(s_2)$. At $t = 0$, we have $h(s,t)=1$ for all $s$. Thus, for $f(s_1)>f(s_2)$, we have

(4.29)\begin{equation} h^{PA}(s_1,t) < h^{PA}(s_2,t),\quad \forall t > 0 \end{equation}

and an analogous argument using (4.27) gives that

(4.30)\begin{equation} h^{S}(s_1,t) < h^{S}(s_2,t),\quad \forall t > 0. \end{equation}

Then, we deduce that the location of minimum thickness for both $h^{PA}$ and $h^{S}$ at any time $t>0$ is given by the Lagrangian point $s_m$ with maximum value for $f$ (explicitly, $s_m:=\arg \max _{s_0} f(s_0)$). In particular, the point of minimum thickness for the purely asymmetric case (PA) and symmetric case (S) will both be at the Lagrangian point $s_0 = s_m$. It should be noted that although the rupture point of the purely asymmetric case (PA) and symmetric case (S) is the same in Lagrangian coordinates, they need not be the same in Eulerian coordinates.

4.6.2. Rate of thinning

Another natural question to ask would be whether the purely asymmetric case (PA) thins faster/slower and, hence, pinches off earlier/later than the symmetric case (S). In order to analyse this problem, we consider two cases where (i) $\sigma = \sigma (\varGamma )$ is convex and (ii) $\sigma = \sigma (\varGamma )$ is concave (see figure 6). From (4.26) and (4.27), we have

(4.31)\begin{equation} \left.\frac{{\rm D}\left(h^{PA}-h^S\right)}{{\rm D}t}\right|_{(s,t)} = \frac{1}{4}\left(\sigma(f(s)h^{PA}(s,t)) - 2\sigma\left(\frac{1}{2}f(s)h^S(s,t)\right) \right),\end{equation}

which we now analyse. From the definition of convexity and concavity, if $\sigma = \sigma (\varGamma )$ is convex, then for $a \in \mathbb {R}$,

(4.32)\begin{equation} \frac{\sigma(a)+\sigma(0)}{2}\geq \sigma\left(\frac{1}{2}a\right) \end{equation}

and if $\sigma = \sigma (\varGamma )$ is concave, then

(4.33)\begin{equation} \frac{\sigma(a)+\sigma(0)}{2}\leq \sigma\left(\frac{1}{2}a\right). \end{equation}

Furthermore, since $\sigma (\varGamma )$ is monotonically decreasing and $\sigma (0)=0$, we have that if $\sigma = \sigma (\varGamma )$ is convex and $h^{PA}(s,t)\geq h^S(s,t)$, then

(4.34)\begin{equation} \sigma(f(s)h^{PA}(s,t)) \geq 2\sigma\left(\frac{1}{2}f(s)h^{PA}(s,t)\right) \geq 2\sigma\left(\frac{1}{2}f(s)h^{S}(s,t)\right) \end{equation}

and if $\sigma = \sigma (\varGamma )$ is concave and $h^{PA}(s,t) \leq h^S(s,t)$, then

(4.35)\begin{equation} \sigma(f(s)h^{PA}(s,t)) \leq 2\sigma\left(\frac{1}{2}f(s)h^{PA}(s,t)\right) \leq 2\sigma\left(\frac{1}{2}f(s)h^{S}(s,t)\right).\end{equation}

At that initial time, $h^{PA}(s,0)=h^S(s,0)=1$. Thus, from (4.31), (4.34) and (4.35), we deduce that

(4.36)\begin{equation} h^{PA}(s,t)\geq h^S(s,t)\quad \forall t \text{ if } \sigma = \sigma(\varGamma) \text{ is convex } \end{equation}

and

(4.37)\begin{equation} h^{PA}(s,t) \leq h^S(s,t)\quad \forall t \text{ if } \sigma = \sigma(\varGamma) \text{ is concave},\end{equation}

which, in other words, says that the purely asymmetric case (PA) thins slower/faster than the symmetric case (S) if $\sigma = \sigma (\varGamma )$ is convex/concave. The conclusion also holds when considering a general asymmetric case (see Appendix E). The conclusion makes sense physically, because a convex/concave $\sigma = \sigma (\varGamma )$ means that the decrease of surface tension is smaller/greater with a larger concentration of surfactant. The purely asymmetric case (PA) would correspondingly have smaller/larger total Marangoni stress than the symmetric case (S).

In figure 7, for $C = 0.5$. we compare the symmetric case (S) (black curves) and the purely asymmetric case (PA) (red curves) with results reported for (a) $\sigma (\varGamma ) = \varGamma (\varGamma -2)$ for $\varGamma \in [0,1]$ at $t = 10$ and (b) $\sigma (\varGamma )= \log (1-{\varGamma }/{2})$ at $t = 20$. In both cases, the initial surfactant concentration is given by $\varGamma _{I+}(x) = \textrm {e}^{-x^2}$ and $\varGamma _{I-}(x) = 0$. The curves plotted are the sheet profile $H\pm \tfrac {1}{2}h$. The isotherm $\sigma = \varGamma (\varGamma -2)$ is chosen to represent some convex isotherm (strictly monotonically decreasing for $\varGamma \in [0,1]$). The isotherm $\sigma (\varGamma )= \log (1-{\varGamma }/{2})$ is chosen to represent an example of a concave isotherm with motivation from the Langmuir isotherm. As expected from (4.36), we see that in figure 7(a) the symmetric case pinches faster than the purely asymmetric case. Similarly, as expected from (4.37), we see in figure 7(b) that the symmetric case pinches more slowly than the asymmetric case. A small point to note is that figure 7 is given in Eulerian coordinates (whereas (4.36) and (4.37) are results in Lagrangian coordinates).

Figure 7. The symmetric case (S) (black curves) and the purely asymmetric case (PA) (red curves) are plotted for (a) $\sigma (\varGamma ) = \varGamma (\varGamma -2)$ at $t = 10$ and (b) $\sigma (\varGamma )= \log (1-{\varGamma }/{2})$ at $t = 20$. The curves plotted are the sheet profile $H\pm \tfrac {1}{2}h$. The spatial variable $x$ is in Eulerian coordinates. In both cases, the initial surfactant concentration is given by $\varGamma _{I+}(x) = \textrm {e}^{-x^2}$ and $\varGamma _{I-}(x) = 0$. $C = 0.5$. For the non-dimensionalisation, see (4.1).

4.7. Early time behaviour

In this subsection, we discuss the early time behaviour for the initial conditions (4.5) and boundary conditions (4.6). It is common in low-Reynolds-number flow analyses to consider inertia at early times in order to analyse how the quasi-static (inertialess) balance is obtained. For the case of liquid sheets and threads, see Buckmaster, Nachman & Ting (Reference Buckmaster, Nachman and Ting1975) and Howell (Reference Howell1996). For many problems, the exact details of the shorter timescale is not as interesting as the later time thinning and pinching as described in §§ 4.1–4.6. Thus, we keep discussions in this section brief. In particular, we only consider $\epsilon ^2 \ll Re \ll 1$. For even smaller $ Re $, analogous approaches can be considered, but the exact equations are different (see Howell Reference Howell1996). In this subsection, we work in Eulerian coordinates.

The main timescale (dimensional) for thinning as described in prior sections is $t_{qs}= {\epsilon \mu \mathcal {L}}/{\Delta \varSigma }$. There are two shorter timescales. The first is the timescale (dimensional) $ Re t_{qs}$, which evolves the horizontal velocity $\bar {u}$, and the second is the timescale (dimensional) $\sqrt { Re }t_{qs}$, which evolves the centreline $H$. Note that since $ Re \ll 1$, we have separation of three timescales $ Re t_{qs} \ll \sqrt { Re }t_{qs} \ll t_{qs}$, which is a familiar concept in asymptotics (Hinch Reference Hinch1991). It will be shown in §§ 4.7.1 and 4.7.2 that on the timescale $ Re t_{qs}$, the horizontal velocity $\bar {u}$ evolves as a diffusion equation and on the timescale $\sqrt { Re } t_{qs}$ the centreline $H$ evolves as a wave equation. For a summary of the timescales, see figure 8. The evolution of $\bar {u}$ on the timescale $ Re t_{qs}$ is an inner solution for the centreline evolution $H$ on the timescale $\sqrt { Re } t_{qs}$, which, in turn, is an inner solution for the quasi-static evolution as considered in §§ 4.1–4.6.

Figure 8. Summary of timescales in § 4. For $t \ll \textit {O}( Re t_{qs})$, the flow is not extensional (red) (see Appendix B). For $t = \textit {O}( Re t_{qs})$ or greater, the flow is extensional (green) (see Appendix B). Three timescales are discussed in § 4: $t = \textit {O}( Re t_{qs})$ (§ 4.7.1), $t = \textit {O}(\sqrt { Re }t_{qs})$ (§ 4.7.2) and $t = \textit {O}(t_{qs})$ (§§ 4.1–4.6).

4.7.1. Diffusion equation for horizontal velocity $\bar {u}$

On the dimensional timescale $ Re t_{qs}$, the non-dimensionalisation is now given by

(4.38)\begin{align} \left. \begin{array}{@{}c@{}} x = \mathcal{L} \tilde{x}, \quad u = \dfrac{\Delta \varSigma}{\epsilon \mu} \tilde{u}, \quad H = \epsilon \mathcal{L} \tilde{H}, \quad t = Re \dfrac{\epsilon \mu\mathcal{L}}{\Delta \varSigma}\tilde{T}_1, \quad y = \epsilon \mathcal{L} \tilde{y}, \quad w = \dfrac{1}{ Re }\dfrac{\Delta \varSigma}{\mu} \tilde{W}_1,\\ h = \epsilon \mathcal{L} \tilde{h}, \quad p = \dfrac{\Delta \varSigma}{\epsilon \mathcal{L}} \tilde{p} ,\quad \sigma = \varSigma +\Delta\varSigma \tilde{\sigma}, \end{array} \right\} \end{align}

and we may apply the same expansion techniques as shown in § 3. Then (again dropping the tildes), the analogues to (3.25), (3.26), (3.27), (3.28) and (4.4) are given, respectively, by

(4.39)$$\begin{gather} \frac{\partial h}{\partial T_1} = 0, \end{gather}$$

(4.40)$$\begin{gather}\frac{\partial H}{\partial T_1} = \bar{W}_1, \end{gather}$$

(4.41)$$\begin{gather}h\frac{\partial \bar{u}}{\partial T_1} = \left(\frac{\partial \sigma_+}{\partial x} + \frac{\partial \sigma_-}{\partial x}\right) + 4 \frac{\partial}{\partial x}\left(h\frac{\partial \bar{u}}{\partial x} \right), \end{gather}$$

(4.42)$$\begin{gather}h \frac{\partial \bar{W}_1}{\partial T_1} =0, \end{gather}$$

(4.43)$$\begin{gather}\frac{\partial \varGamma_{{\pm}}}{\partial T_1} = 0. \end{gather}$$

We can immediately deduce from the initial conditions (4.5) that

(4.44)\begin{equation} h = 1, H = 0, \bar{W}_1 = 0,\varGamma_{{\pm}} = \varGamma_{I\pm}, \end{equation}

which then gives that

(4.45)\begin{equation} \frac{\partial \bar{u}}{\partial T_1} - 4 \frac{\partial^2 \bar{u}}{\partial x^2}= \frac{\partial}{\partial x}\left(\sigma(\varGamma_{I+})+\sigma(\varGamma_{I-})\right).\end{equation}

Thus, we see that the evolution equation for $\bar {u}$ is a diffusion equation with diffusion coefficient $4$ and forcing given by the right-hand side of (4.45). From standard considerations of the diffusion equation (see Appendix F), it can then be deduced that

(4.46)\begin{equation} \lim_{T_1 \rightarrow \infty}\frac{\partial \bar{u}}{\partial x} ={-} \frac{1}{4}\left(\sigma(\varGamma_{I+})+\sigma(\varGamma_{I-})\right), \end{equation}

which satisfies the quasi-static balance (4.2) as expected. A subtle point to be made here is that the initial condition (4.5) is compatible with the boundary condition (4.6) because the early time dynamics discussed here respects the far-field condition ${\partial \bar {u}}/{\partial x}= 0$ (explicitly, $\lim _{x \rightarrow \pm \infty } \lim _{T_1 \rightarrow \infty }({\partial \bar {u}}/{\partial x}) = 0$).

4.7.2. Wave equation for centreline H

On the dimensional timescale $\sqrt { Re } t_{qs}$, the non-dimensionalisation is now given by

(4.47)\begin{equation} \left. \begin{array}{@{}c@{}} x = \mathcal{L} \tilde{x}, \quad u = \dfrac{\Delta \varSigma}{\epsilon \mu} \tilde{u}, \quad H = \epsilon \mathcal{L} \tilde{H}, \quad t = \sqrt{ Re } \dfrac{\epsilon \mu \mathcal{L}}{\Delta \varSigma}\tilde{T}_2, \quad y = \epsilon \mathcal{L} \tilde{y},\\ w = \dfrac{1}{\sqrt{ Re }}\dfrac{\Delta \varSigma}{\mu} \tilde{W}_2,\quad h = \epsilon \mathcal{L} \tilde{h}, \quad p = \dfrac{\Delta \varSigma}{\epsilon \mathcal{L}}\tilde{p} ,\quad \sigma = \varSigma +\Delta\varSigma \tilde{\sigma}, \end{array} \right\} \end{equation}

and we may apply the same expansion techniques as shown in § 3. Then, the analogues to (3.25), (3.26), (3.27), (3.28) and (4.4) are given, respectively, by

(4.48)$$\begin{gather} \frac{\partial h}{\partial T_2} = 0, \end{gather}$$

(4.49)$$\begin{gather}\frac{\partial H}{\partial T_2} = \bar{W}_2, \end{gather}$$

(4.50)$$\begin{gather}0 = \left(\frac{\partial \sigma_+}{\partial x} + \frac{\partial \sigma_-}{\partial x}\right) + 4 \frac{\partial}{\partial x}\left(h\frac{\partial \bar{u}}{\partial x} \right), \end{gather}$$

(4.51) $$\begin{gather} \frac{\partial \bar{W}_2}{\partial T_2} =\frac{1}{2}\frac{\partial^2 }{\partial x^2}\left(h(\sigma_{+}-\sigma_-)\right)+\frac{\partial}{\partial x}\left(\frac{\partial H}{\partial x}\left(\frac{2}{C}+\left(\sigma_{+}+\sigma_-\right)\right)\right) \nonumber\\ +\ 4 \frac{\partial}{\partial x}\left(h \frac{\partial H}{\partial x}\frac{\partial \bar{u}}{\partial x}\right), \end{gather}$$

(4.52)$$\begin{gather}\frac{\partial \varGamma_{{\pm}}}{\partial T_2} = 0. \end{gather}$$

Regarding the shorter timescale $T_1$ (in dimensional variables, $ Re t_{qs}$) of the evolution for $\bar {u}$ described in § 4.7.1 as the inner solution, we deduce from (4.46) that on the timescale $T_2$ of centreline evolution,

(4.53)\begin{equation} \frac{\partial \bar{u}}{\partial x} ={-} \frac{1}{4}\left(\sigma(\varGamma_{I+})+\sigma(\varGamma_{I-})\right), \end{equation}

along with

(4.54)\begin{equation} h = 1\quad\text{and}\quad \varGamma_{{\pm}} = \varGamma_{I\pm}, \end{equation}

which can be substituted into (4.49) and (4.51) to get

(4.55)\begin{equation} \frac{\partial^2 H}{\partial T_2^2} - \frac{2}{C}\frac{\partial^2 H}{\partial x^2} = \frac{1}{2}\frac{\partial^2}{\partial x^2}\left(\sigma(\varGamma_{I+})-\sigma(\varGamma_{I-})\right). \end{equation}

Equation (4.55) is a wave equation with wavespeed $c = \sqrt {{2}/{C}}$ and forcing term independent of $T_2$ given by the right-hand side of (4.55). In dimensional variables, the wavespeed is $({\Delta \varSigma }/{\epsilon \mu }) Re ^{-{1}/{2}}\sqrt {{2}/{C}}$. Thus, the solution can be written as

(4.56)\begin{equation} H ={-} \frac{1}{4c^2}\left(2g(x)-g(x-cT_2)-g(x+cT_2)\right),\end{equation}

where $g(x):=\sigma (\varGamma _{I+}(x))-\sigma (\varGamma _{I-}(x))$. Equation (4.56) is the famous d'Alembert's solution for the wave equation and can be interpreted as having two travelling waves, one travelling to the left with wavespeed $c$ and one travelling to the right with wavespeed $c$. In particular, we may deduce that

(4.57)\begin{equation} \lim_{T_2 \rightarrow \infty} H ={-} \frac{C}{4}\left(\sigma(\varGamma_{I+})- \sigma(\varGamma_{I-})\right), \end{equation}

which satisfies the quasi-static balance (4.3) as expected. Again, it should be noted that the initial condition (4.5) is, thus, compatible with the boundary condition (4.6) since the early time dynamics discussed here respects the far field condition $H = 0$.

5.1. Non-dimensional equations

Non-dimensionalise the equations according to (3.1) with $\mathcal {U}=\mathcal {L}\alpha$ (in other words, consider the timescale of the background exponentially thinning flow), where $\epsilon = {h_I}/{\mathcal {L}}$. In addition, non-dimensionalise the surfactant concentration with a characteristic surfactant concentration $\varGamma _c$. Explicitly,

(5.1)\begin{equation} \left. \begin{array}{@{}c@{}} x = \mathcal{L} \tilde{x}, \quad u = \mathcal{L} \alpha \tilde{u}, \quad H = \epsilon \mathcal{L} \tilde{H}, \quad t = \dfrac{1}{\alpha}\tilde{t}, \quad y = \epsilon \mathcal{L} \tilde{y}, \quad w = \epsilon \alpha \mathcal{L} \tilde{w}, \\ h = \epsilon \mathcal{L} \tilde{h}, \quad \sigma = \varSigma+ \Delta \varSigma \tilde{\sigma}, \quad \varGamma = \varGamma_c \tilde{\varGamma}, \end{array} \right\} \end{equation}

We then have parameters $\textit {M}={\Delta \varSigma }/{\epsilon \mu \mathcal {L}\alpha }$, $ Re ={\rho \mathcal {L}^2\alpha }/{\mu }$ and $C={\epsilon \mu \mathcal {L} \alpha }/{\varSigma }$. Again, the tilde $\widetilde {(\ )}$ is dropped. We consider $M = \textit {O}(1)$ or less and $C = \textit {O}(1)$ or more so that we have extensional flow (see § 3.2). The horizontal length scale $\mathcal {L}$ is once again given by the width of variation of the initial surfactant concentration.

Equations for $h$ and $H$ are given by (3.25) and (3.26). The equations for $\bar {u}$ and $\bar {w}$ are given by (3.27) and (3.28) with the inertia term neglected (note that unlike § 4, the constant $M$ is not necessarily $1$). Diffusion is ignored once again by assuming a large enough Péclet number ${\mathcal {L}^2 \alpha }/{D }$ (diffusion coefficient $D$). Thus, $\varGamma$ is given by (4.4).

5.2. Initial and boundary conditions

For the initial conditions, we consider a uniform sheet with constant extension, along with some initial surfactant distribution $\varGamma _{I\pm }(x)$ where $\varGamma _{\pm } \rightarrow 0$ as $x \rightarrow \pm \infty$. Far-field conditions are taken for the boundary condition. See figure 9 for a particular example where $\varGamma _{I+}=\textrm {e}^{-x^2}$ and $\varGamma _{I-} = 0$, which is considered in § 5.4. Explicitly,

(5.2)\begin{equation} h = 1,\quad H = 0,\quad \bar{u} = x,\quad \bar{w} = 0,\quad \varGamma_{{\pm}} = \varGamma_{I\pm} \quad \text{at } t = 0\end{equation}

and

(5.3)\begin{equation} h = {\rm e}^{{-}t},\quad H = 0,\quad \frac{\partial \bar{u}}{\partial x}= 1,\quad \bar{w} = 0, \quad \varGamma_{{\pm}} = 0\quad \text{at } x ={\pm} \infty, \end{equation}

where it should be noted here that $x = \pm \infty$ more rigorously means some $1 \ll x \ll Re ^{-1}$ (in order for the inertia terms to be neglected, $ Re \bar {u} ({\partial \bar {u}}/{\partial x}) = Re x \ll 1$).

Figure 9. Initial set-up considered in § 5.4. An initial Gaussian surfactant distribution $\varGamma _{I+}(x) = \textrm {e}^{-x^2}$, $\varGamma _{I-}(x) = 0$ is deposited on top of a uniform sheet $-\tfrac {1}{2} \leq z \leq \tfrac {1}{2}$ at $t = 0$ with initial condition for horizontal velocity given by $\bar {u}_I(x) = x$. The horizontal direction is given by the $x$ axis and the vertical direction is given by the $z$ axis. The red arrows denote the direction of horizontal flow. See (5.1) for the non-dimensionalisation.

5.3. Lagrangian coordinates

We may proceed with the same technique as in § 4.3 to deduce the equivalent equations to (4.7), (4.8), (4.12), (4.13), (4.14) and (4.15):

(5.4)$$\begin{gather} \frac{\partial \bar{u}}{\partial x} ={-}\frac{M}{4h(x,t)} \left(\sigma(\varGamma_+(x,t))+\sigma(\varGamma_-(x,t))\right) +\frac{{\rm e}^{{-}t}}{h(x,t)}, \end{gather}$$

(5.5)$$\begin{gather}H(x,t) ={-}\frac{M}{2}\frac{h(\sigma(\varGamma_+(x,t))-\sigma(\varGamma_-(x,t)))}{\left(\dfrac{2}{C}+4{\rm e}^{{-}t}\right)}, \end{gather}$$

(5.6)$$\begin{gather}\varGamma_{{\pm}}(s,t) = \varGamma_{I\pm}(s)h(s,t), \end{gather}$$

(5.7)$$\begin{gather}\left.\frac{{\rm D}h}{{\rm D}t}\right|_{(s,t)} =\frac{M}{4}\left(\sigma(\varGamma_{I+}(s)h(s,t))+\sigma(\varGamma_{I-}(s)h(s,t))\right) - {\rm e}^{{-}t}, \end{gather}$$

(5.8)$$\begin{gather}\bar{w}(s,t) = \left.\frac{{\rm D}H}{{\rm D}t}\right|_{(s,t)} , \end{gather}$$

(5.9)$$\begin{gather}\left.\frac{\partial x}{\partial s}\right|_{(s,t)} = \frac{1}{h(s,t)}. \end{gather}$$

5.4. Example: linear isotherm, Gaussian deposition on top

In this subsection, we consider the example where the initial surfactant distribution given by $\varGamma _{I+}(s) = \textrm {e}^{-s^2}$ and $\varGamma _{I-}(s) = 0$ and the surface tension isotherm is linear $\sigma (\varGamma ) = - \varGamma$ (see figure 9 for the set-up). Then, we have

(5.10)\begin{equation} \left.\frac{{\rm D}h}{{\rm D}t}\right|_{(s,t)} ={-} K(s) h(s,t) - {\rm e}^{{-}t}, \end{equation}

where $K(s):= ({M}/{4})\textrm {e}^{-s^2}$. Using the initial condition $h(s,0) = 1$, we have the analytic solution

(5.11)\begin{equation} h(s,t) = \left\{ \begin{array}{@{}ll} \dfrac{{\rm e}^{{-}t}}{1-K(s)}-\dfrac{K(s)}{1-K(s)}{\rm e}^{{-}K(s)t}, & K(s) \neq 1, \\ (1-t){\rm e}^{{-}t}, & K(s) = 1, \end{array} \right. \end{equation}

where we note that the point $K(s)=1$ is not a discontinuity since

(5.12)\begin{equation} \lim_{k \rightarrow 1} \left( \frac{{\rm e}^{{-}t}}{1-k}-\frac{k}{1-k}{\rm e}^{{-}kt}\right) = (1-t){\rm e}^{{-}t}. \end{equation}

Noting that $x(0,t)=0$ by symmetry, we have from (5.9) that

(5.13)\begin{equation} x(s,t) = \int_0^{s}\frac{1}{h(s',t)} \,{\rm d}s'.\end{equation}

The equations for $H$, $\bar {u}$, $\bar {w}$ and $\varGamma _{\pm }$ are given by

(5.14)\begin{gather} H(s,t) = \frac{2K(s) h(s,t)^2}{\left(\dfrac{2}{C}+ 4 {\rm e}^{{-}t}\right)}, \end{gather}

(5.15)\begin{gather} \bar{u}(s,t) = \int_0^{s} \frac{K(s')}{h(s',t)} + \frac{{\rm e}^{{-}t}}{h^2(s',t)} \,{\rm d}s', \end{gather}

(5.16)\begin{gather} \bar{w}(s,t) ={-}\frac{4K(s) h(s,t)}{\left(\dfrac{2}{C}+ 4 {\rm e}^{{-}t}\right)}\left(K(s )h(s,t) +{\rm e}^{{-}t}\right) + \frac{8{\rm e}^{{-}t}K(s) h(s,t)^2}{\left(\dfrac{2}{C}+ 4 {\rm e}^{{-}t}\right)^2}, \end{gather}

(5.17)\begin{gather} \varGamma_+(s,t) = \frac{4}{M}K(s)h(s,t),\quad \varGamma_-(s,t) =0. \end{gather}

5.4.1. Description of the solution

First, we note from (5.11) that $h(s,t)=0$ when $t = ({1}/({K(s)-1}))\log (K(s))$ for $K(s)\neq 1$ or when $t = 1$ for $K(s)=1$. Thus, rupture occurs at time $t = t_r$ and location $s = s_r$ where

(5.18)\begin{equation} t_r = \min_{s} \left\{ \begin{array}{@{}ll@{}} \dfrac{1}{K(s)-1}\log(K(s)), & K(s) \neq 1, \\ 1, & K(s) = 1,\end{array} \right. \end{equation}

and

(5.19)\begin{equation} s_r = \arg \min_{s} \left\{ \begin{array}{@{}ll@{}} \dfrac{1}{K(s)-1}\log(K(s)), & K(s) \neq 1, \\ 1, & K(s) = 1, \end{array} \right. \end{equation}

which, by noting that the function $f(k) :=\log (k)/(k-1)$ is a strictly decreasing function of $x>0$, can be evaluated as

(5.20)\begin{equation} t_r = \left\{ \begin{array}{@{}ll@{}} \dfrac{1}{\dfrac{M}{4}-1}\log\left(\dfrac{M}{4}\right), & M \neq 4, \\ 1, & M = 4, \end{array} \right. \end{equation}

and

(5.21)\begin{equation} s_r = 0. \end{equation}

There is no discontinuity once again since $\lim _{M \rightarrow 4} ({1}/({{M}/{4}-1}))\log ({M}/{4}) = 1$.

It is interesting to note that for any $M>0$, finite time rupture now occurs (in contrast to the case without surfactants where the sheet is exponentially thinning $h = \textrm {e}^{-t}$). The fact that $s_r = 0$ follows once again from the more general notion that rupture occurs at the Lagrangian point of maximum $\varGamma$.

Previously, we analysed where and when the sheet ruptures. We may also analyse how the sheet ruptures. Again, a point of caution is that this paper does not include van der Waals forces. As shown in Appendix G, we have for $s \ll 1$ and $0< t^*:=t_r-t \ll 1$ that

(5.22)\begin{equation} h(s,t^*) = c_1 t^* + c_2 s^2\end{equation}

and

(5.23)\begin{equation} x(s,t^*) = \frac{1}{\sqrt{c_1c_2 t^*}}\arctan\left(\frac{s}{\sqrt{\dfrac{c_1}{c_2}t^*}}\right), \end{equation}

where $c_1$, $c_2$ are explicit functions of $M$, given by (G2). Equations (5.22) and (5.23) can be rewritten as self-similar solutions

(5.24)\begin{equation} h(\eta,t^*) = c_1t^*\left(1+\eta^2\right) \end{equation}

and

(5.25)\begin{equation} x(\eta,t^*)= \frac{1}{\sqrt{c_1c_2t^*}}\arctan(\eta), \end{equation}

where $\eta :=s(({c_1}/{c_2})t^*)^{-{1}/{2}}$ is the self-similarity variable. Inverting (5.23) and substituting into (5.22) gives that near the pinch ($s \ll 1$, $t^* \ll 1$), $h$ is given in Eulerian coordinates by

(5.26)\begin{equation} h(x,t^*) = c_1 t^*\left(1+ \tan^2 \left(x\sqrt{c_1c_2t^*}\right)\right). \end{equation}

An example with $M = 4$ and $C = 0.5$ (representative of parameter values with $M = \textit {O}(1)$ and $C =\textit {O}(1)$) is given in figure 10. The value $M =4$ has the nice property from (5.20) that the sheet ruptures at $t = 1$. Figure 10(a) shows the bottom of the sheet $H-\tfrac {1}{2}h$, the centreline $H$, the top of the sheet $H+ \tfrac {1}{2}h$, and the surfactant distribution $\varGamma _{+}$ at the initial time. Again, by initial time, we mean $t = 0$ where a quasi-static balance has been achieved and so, in particular, the centreline is bent (see §§ 5.1 and 5.5). Figure 10(b) shows the velocity profiles $\bar {u}-x$ and $\bar {w}$ at the initial time. The value $\bar {u}-x$ is plotted rather than $\bar {u}$ for better visualisation. Similarly, figure 10(c,d) report the variables at $t = 0.5$ and figure 10(e, f) report the variables at $t = 0.9$ (which is close to the rupture time $t_r = 1$). The times $t = 0.5$ and $0.9$ were chosen to represent mid and late times, respectively, for the thinning problem.

Figure 10. Time evolution with an imposed extensional flow. Plots of (5.11), (5.14), (5.15), (5.16), (5.17) and (5.13) for $M=4$ and $C=0.5$. The spatial variable $x$ is in Eulerian coordinates. The three curves from bottom to top in (a) are $H - \tfrac {1}{2}h$ (black), $H$ (blue), $H + \tfrac {1}{2}h$ (coloured with surfactant concentration $\varGamma$) at the initial time and (b) shows the velocity profiled $\bar {u}-x$ (blue) and $\bar {w}$ (red) at the initial time. Similarly, (c,d) report the variables at $t = 0.5$ and (e, f) report the variables at $t = 0.9$. For the non-dimensionalisation, see (4.1).

The non-uniform surfactant distribution means that there are Marangoni effects, where surfactants spread from a region of high concentration to low concentration (see figure 10a,c,e), since low-concentration regions have higher surface tension. As result of the horizontal velocity $\bar {u}$ induced by the spreading of surfactants, $h$ decreases and the sheet thins, in addition to the background exponentially thinning given by $\textrm {e}^{-t}$; see figure 10(a,c,e).

5.5. Early time behaviour

Again, we only consider $\epsilon ^2 \ll Re \ll 1$. There are three timescales once again, given by $ Re \alpha ^{-1} \ll \sqrt { Re }\alpha ^{-1} \ll \alpha ^{-1}$. On the timescale $ Re \alpha ^{-1}$, the horizontal velocity $\bar {u}$ evolves as a diffusion equation. Proceeding as in § 4.7.1, we deduce that

(5.27)\begin{equation} \frac{\partial \bar{u}}{\partial T_1} - 4 \frac{\partial^2 \bar{u}}{\partial x^2}= M \frac{\partial}{\partial x}\left(\sigma(\varGamma_{I+})+\sigma(\varGamma_{I-})\right)\end{equation}

with the limit $T_1 \rightarrow \infty$ given by

(5.28)\begin{equation} \lim_{T_1 \rightarrow \infty}\frac{\partial \bar{u}}{\partial x} = 1 - \frac{M}{4}\left(\sigma(\varGamma_{I+})+\sigma(\varGamma_{I-})\right). \end{equation}

On the timescale $\sqrt { Re }\alpha ^{-1}$ the centreline $H$ evolves as a wave equation. Proceeding as in § 4.7.2, we deduce that

(5.29)\begin{equation} \frac{\partial^2 H}{\partial T_2^2} - \left(4 + \frac{2}{C}\right)\frac{\partial^2 H}{\partial x^2} = \frac{M}{2}\frac{\partial^2}{\partial x^2}\left(\sigma(\varGamma_{I+})-\sigma(\varGamma_{I-})\right),\end{equation}

which is the wave equation with wavespeed $c = \sqrt {(4 + {2}/{C})}$ and forcing term independent of $T_2$ given by the right-hand side of (5.29). In dimensional variables, the wavespeed is $\alpha \mathcal {L} Re ^{-{1}/{2}}\sqrt {(4 + {2}/{C})}$. Thus, the solution can be written as

(5.30)\begin{equation} H ={-} \frac{M}{4c^2}\left(2g(x)-g(x-cT_2)-g(x+cT_2)\right), \end{equation}

where $g(x):=\sigma (\varGamma _{I+}(x))-\sigma (\varGamma _{I-}(x))$. In particular, we may deduce that

(5.31)\begin{equation} \lim_{T_2 \rightarrow \infty} H ={-} \frac{M}{2\left(4 + \dfrac{2}{C}\right)}\left(\sigma(\varGamma_{I+})- \sigma(\varGamma_{I-})\right). \end{equation}