3.1 Solution for a single mode

To this point, the derivation has considered an arbitrary initial surfactant distribution $\tilde{\unicode[STIX]{x1D6E4}}_{0}(\tilde{x})$ , and we now restrict our attention to a particular distribution: a single sinusoidal mode. Though a pure sinusoid is an unlikely initial condition for a physical system, the solution for any periodic initial condition can be expressed as a series of the solutions for individual modes due to the linearity of the governing system of equations. Furthermore, the solution for a pure sinusoid is useful because the fundamental mode decays slowest, as we will show, and therefore determines the late-time dynamics. Choosing the length scale $\unicode[STIX]{x1D706}$ to be the wavelength of the initial sinusoidal variation in the surfactant concentration, the initial distribution $\tilde{\unicode[STIX]{x1D6E4}}_{0}(\tilde{x})$ is

(3.11) $$\begin{eqnarray}\tilde{\unicode[STIX]{x1D6E4}}_{0}(\tilde{x})=\sin (2\unicode[STIX]{x03C0}\tilde{x}).\end{eqnarray}$$

We now take the Fourier and Laplace transforms of the surfactant evolution equation (3.1a ) to obtain

(3.12) $$\begin{eqnarray}s\breve{\unicode[STIX]{x1D6E4}}(k,s)-\hat{\unicode[STIX]{x1D6E4}}_{0}(k)+\text{i}k\breve{u}_{s}(k,s)=0,\end{eqnarray}$$

where $\breve{u}_{s}(k,s)=\breve{u}(k,0,s)$ and $\hat{\unicode[STIX]{x1D6E4}}_{0}(k)=\hat{\unicode[STIX]{x1D6E4}}(k,0)$ . Substituting (3.9a ) and (3.10), we arrive at the Laplace-domain solution for the surfactant evolution:

(3.13) $$\begin{eqnarray}\breve{\unicode[STIX]{x1D6E4}}(k,s)=\frac{1}{\displaystyle s+\frac{|k|^{3}}{\mathit{Re}\,s}\left(\sqrt{1+\frac{\mathit{ Re }\,s}{k^{2}}}-1\right)}\hat{\unicode[STIX]{x1D6E4}}_{0}(k).\end{eqnarray}$$

For the initial condition we are interested in, a single sinusoidal mode, we have $\hat{\unicode[STIX]{x1D6E4}}_{0}(k)=-\text{i}\unicode[STIX]{x03C0}\,\text{sgn}(k)\,\unicode[STIX]{x1D6FF}(|k|-2\unicode[STIX]{x03C0})$ , where $\unicode[STIX]{x1D6FF}(\cdot )$ is the Dirac delta function and $\text{sgn}(\cdot )$ is the sign function, and therefore need to only consider $|k|=2\unicode[STIX]{x03C0}$ . The resulting solution is

(3.14) $$\begin{eqnarray}\breve{\unicode[STIX]{x1D6E4}}(k,s)=\frac{-\text{i}\unicode[STIX]{x03C0}\,\text{sgn}(k)\,\unicode[STIX]{x1D6FF}(|k|-2\unicode[STIX]{x03C0})}{\displaystyle s+\frac{8\unicode[STIX]{x03C0}^{3}}{\mathit{Re}\,s}\left(\sqrt{1+\frac{\mathit{ Re }\,s}{4\unicode[STIX]{x03C0}^{2}}}-1\right)}.\end{eqnarray}$$

To gain first a qualitative understanding of the surfactant dynamics, we examine the poles of (3.14). As illustrated in figure 2, the dynamics can be classified into one of five cases depending on the value of the Reynolds number:

Figure 2.

Poles of the Laplace-domain solution for the surfactant evolution (3.14). The blue curve shows the real part of the poles, and the red curve shows their imaginary part. The dotted lines show the components of the poles in the limit of large $\mathit{Re}$ (3.17).

1. (i) $\mathit{Re}\rightarrow 0$ . In this limit, which describes dominantly viscous systems, the pole is $s=-\unicode[STIX]{x03C0}$ . In the absence of the term on the left-hand side of (2.10), we have a Stokes flow. The complete solution for this case can be derived either directly from the governing equations or using the $\mathit{Re}\rightarrow 0$ limit of the solutions that have already been derived, which we present in § 3.2. The Stokes flow solution for a two-dimensional surface and a finite fluid depth is also available (Masoud & Shelley Reference Masoud and Shelley2014).

2. (ii) $\mathit{Re}/(2\unicode[STIX]{x03C0})<1$ , for which only one negative real pole exists.

3. (iii) $1<\mathit{Re}/(2\unicode[STIX]{x03C0})<32/27$ , for which two negative real poles exist.

4. (iv) $\mathit{Re}/(2\unicode[STIX]{x03C0})>32/27$ , for which two complex conjugate poles with negative real parts exist.

5. (v) $\mathit{Re}\gg 1$ . In this limit, equation (3.14) simplifies to

   (3.15) $$\begin{eqnarray}\breve{\unicode[STIX]{x1D6E4}}(k,s)=-\text{i}\unicode[STIX]{x03C0}\,\text{sgn}(k)\,\unicode[STIX]{x1D6FF}(|k|-2\unicode[STIX]{x03C0})\,\breve{\unicode[STIX]{x1D703}}(s),\end{eqnarray}$$
   where
   (3.16) $$\begin{eqnarray}\breve{\unicode[STIX]{x1D703}}(s)=\frac{1}{\displaystyle s+\frac{4\unicode[STIX]{x03C0}^{2}}{\mathit{Re}^{1/2}s^{1/2}}}\end{eqnarray}$$
   and has the two poles
   (3.17) $$\begin{eqnarray}s=2\unicode[STIX]{x03C0}\left(\frac{\mathit{Re}}{2\unicode[STIX]{x03C0}}\right)^{-1/3}\left(-\frac{1}{2}\pm \text{i}\frac{\sqrt{3}}{2}\right).\end{eqnarray}$$
   The poles are complex conjugates of each other and have negative real parts. Both the real and imaginary parts of the two poles decay to zero in proportion to $\mathit{Re}^{-1/3}$ .

We focus now on the limit of large $\mathit{Re}$ and invert the Laplace-domain solution, equation (3.15), back to the time domain. First, we note that

(3.18) $$\begin{eqnarray}\displaystyle \breve{\unicode[STIX]{x1D703}}(s) & = & \displaystyle \frac{1}{\displaystyle s+\frac{4\unicode[STIX]{x03C0}^{2}}{\mathit{Re}^{1/2}\,s^{1/2}}}\nonumber\\ \displaystyle & = & \displaystyle \frac{1}{3s^{1/2}}\left[{\displaystyle \frac{1}{s^{1/2}-b_{0}\frac{\displaystyle (2\unicode[STIX]{x03C0})^{2/3}}{\displaystyle \mathit{Re}^{1/6}}}}+{\displaystyle \frac{1}{s^{1/2}-b_{+}\frac{\displaystyle (2\unicode[STIX]{x03C0})^{2/3}}{\displaystyle \mathit{Re}^{1/6}}}}+{\displaystyle \frac{1}{s^{1/2}-b_{-}\frac{\displaystyle (2\unicode[STIX]{x03C0})^{2/3}}{\displaystyle \mathit{Re}^{1/6}}}}\right],\qquad\end{eqnarray}$$

where $b_{0}=-1$ and $b_{\pm }=(1\pm \text{i}\sqrt{3})/2$ are the three cube roots of $-1$ . Of the three terms on the right, the first does not contribute a pole. The inverse Laplace transform of

(3.19) $$\begin{eqnarray}F(s)=\frac{1}{\sqrt{s}\left(\sqrt{s}+\unicode[STIX]{x1D6FD}\right)},\end{eqnarray}$$

where $\unicode[STIX]{x1D6FD}$ is a complex constant, is (Erdélyi et al. Reference Erdélyi, Magnus, Oberhettinger and Tricomi1954, p. 234, equation (7))

(3.20) $$\begin{eqnarray}\mathscr{L}^{-1}\{F(s)\}=f(t)=\exp (\unicode[STIX]{x1D6FD}^{2}t)\text{erfc}(\unicode[STIX]{x1D6FD}\sqrt{t}),\end{eqnarray}$$

and we can use this to invert every term in (3.18). We thus obtain

(3.21) $$\begin{eqnarray}\hat{\unicode[STIX]{x1D6E4}}(k,\tilde{t})=-\text{i}\unicode[STIX]{x03C0}\,\text{sgn}(k)\,\unicode[STIX]{x1D6FF}(|k|-2\unicode[STIX]{x03C0})\,\unicode[STIX]{x1D703}(\tilde{t})\end{eqnarray}$$

and

(3.22) $$\begin{eqnarray}\tilde{\unicode[STIX]{x1D6E4}}(\tilde{x},\tilde{t})=\sin (2\unicode[STIX]{x03C0}\tilde{x})\,\unicode[STIX]{x1D703}(\tilde{t}),\end{eqnarray}$$

where the time-dependent component, $\unicode[STIX]{x1D703}(\tilde{t})$ , which determines how the amplitude of the surfactant distribution decays over time, is

(3.23) $$\begin{eqnarray}\unicode[STIX]{x1D703}(\tilde{t})=\frac{1}{3}\mathop{\sum }_{n=0,+,-}\exp \left(b_{n}^{2}\frac{(2\unicode[STIX]{x03C0})^{4/3}}{\mathit{Re}^{1/3}}\tilde{t}\right)\text{erfc}\left(-b_{n}\frac{(2\unicode[STIX]{x03C0})^{2/3}}{\mathit{Re}^{1/6}}\sqrt{\tilde{t}}\right).\end{eqnarray}$$

For any complex $z$ , $\exp (\overline{z}^{2})\text{erfc}(\overline{z})=\overline{\exp (z^{2})\text{erfc}(z)}$ , where an overbar denotes the complex conjugate. Since $b_{-}=\overline{b_{+}}$ , the $n=+$ and $n=-$ terms of the sum are complex conjugates of each other, and the sum of the three terms is real.

The decay of the variation in surfactant concentration resembles a damped oscillation, as shown in figure 3, but the amplitude of the variation becomes zero only three times and decays monotonically at late times. Determined numerically, the three roots $\tilde{t}_{i}$ of $\unicode[STIX]{x1D703}(\tilde{t})$ are $((2\unicode[STIX]{x03C0})^{4/3}/\mathit{Re}^{1/3})\tilde{t}_{i}=1.6452$ , $5.7437$ and $8.3765$ . Oscillations occur due to the inertia of the fluid that is set in motion by the surface stress distribution. Initially, regions with high surfactant concentration (lower surface tension) expand, while those with low concentration (higher tension) contract. When the surfactant distribution becomes uniform, first at $\tilde{t}_{1}$ , the fluid continues to expand/contract the surfactant distribution, which now creates regions of locally low/high concentration. The surface stresses then reverse the upper layer of the flow, and the surfactant distribution decays again. The surfactant concentration overshoots the equilibrium two more times before gradually decaying to uniformity. The finite number of oscillations distinguishes this system from a typical damped oscillator.

Figure 3.

The decay over time of the amplitude of a sinusoidal surfactant concentration distribution at high $\mathit{Re}$ (3.23). The early time behaviour (3.26) is shown with a dash-dotted line. The inset emphasizes the finite number of roots and the absence of sustained oscillation, with the late-time monotonic decay (3.24) shown as a dashed line.

Using $\exp (z^{2})\text{erfc}(z)\sim (1/z\sqrt{\unicode[STIX]{x03C0}})(1-1/2z^{2})$ for large $z$ , obtained from the first two terms of the asymptotic expansion of $\text{erfc}$ (Olver et al. Reference Olver, Lozier, Boisvert and Clark2010, equation (7.12.1); DLMF 2015), we find that at late times, deviations in surfactant concentration decay monotonically according to the power law

(3.24) $$\begin{eqnarray}\unicode[STIX]{x1D703}(\tilde{t})\sim -\frac{1}{2\sqrt{\unicode[STIX]{x03C0}}}\left(\frac{(2\unicode[STIX]{x03C0})^{4/3}}{\mathit{Re}^{1/3}}\tilde{t}\right)^{-3/2}.\end{eqnarray}$$

This result is shown as the dashed line in the inset of figure 3. Two terms in the approximation are required because $1/b_{0}+1/b_{+}+1/b_{-}=0$ , which makes the contribution of the first term vanish in the sum in equation (3.23).

The factor on $\tilde{t}$ in the solution for the surfactant decay, equation (3.23), shows that in the large $\mathit{Re}$ limit, the decay time scale is not the advective time scale $t^{\ast }=\unicode[STIX]{x1D706}/u^{\ast }$ . Instead, the surfactant decay time scale $t_{\unicode[STIX]{x1D6E4}}^{\ast }$ is the cube root of the product of the viscous time scale $\unicode[STIX]{x1D70C}\unicode[STIX]{x1D706}^{2}/\unicode[STIX]{x1D707}$ and the square of the advective time scale $t^{\ast }=\unicode[STIX]{x1D706}/u^{\ast }$ :

(3.25) $$\begin{eqnarray}t_{\unicode[STIX]{x1D6E4}}^{\ast }=t^{\ast }\mathit{Re}^{1/3}=\left(\frac{\unicode[STIX]{x1D70C}\unicode[STIX]{x1D706}^{4}\unicode[STIX]{x1D707}}{\unicode[STIX]{x1D6FC}^{2}\bar{\unicode[STIX]{x1D6E4}}^{2}}\right)^{1/3}=\left(t^{\ast 2}\frac{\unicode[STIX]{x1D70C}\unicode[STIX]{x1D706}^{2}}{\unicode[STIX]{x1D707}}\right)^{1/3}.\end{eqnarray}$$

One may still interpret this as an advective time scale with characteristic velocity $u^{\ast }/\mathit{Re}^{1/3}$ instead of $u^{\ast }$ .

The surfactant decay time scale also follows from a simple scaling argument. From the surfactant transport equation (3.1a ), we have $\unicode[STIX]{x2202}\tilde{\unicode[STIX]{x1D6E4}}/\unicode[STIX]{x2202}\tilde{t}\sim \unicode[STIX]{x2202}\tilde{u} _{s}/\unicode[STIX]{x2202}\tilde{x}\sim \tilde{u} _{s}$ , and the surface stress boundary condition (2.8) provides $\unicode[STIX]{x2202}\tilde{u} /\unicode[STIX]{x2202}{\tilde{y}}\sim \unicode[STIX]{x2202}\tilde{\unicode[STIX]{x1D6E4}}/\unicode[STIX]{x2202}\tilde{x}\sim \tilde{\unicode[STIX]{x1D6E4}}$ for ${\tilde{y}}\rightarrow 0^{+}$ . Making the approximation that $\unicode[STIX]{x2202}\tilde{u} /\unicode[STIX]{x2202}{\tilde{y}}\sim \tilde{u} _{s}/\tilde{\unicode[STIX]{x1D6FF}}$ , where $\tilde{\unicode[STIX]{x1D6FF}}\ll 1$ is a boundary layer thickness that grows over time as $\tilde{\unicode[STIX]{x1D6FF}}\sim (\tilde{t}/\mathit{Re})^{1/2}$ (this follows from the momentum equation (3.1b ) and $\mathit{Re}\,\unicode[STIX]{x2202}\tilde{u} /\unicode[STIX]{x2202}\tilde{t}\sim \unicode[STIX]{x2202}^{2}\tilde{u} /\unicode[STIX]{x2202}{\tilde{y}}^{2}$ ), we obtain $\tilde{u} _{s}\sim \tilde{\unicode[STIX]{x1D6E4}}\tilde{\unicode[STIX]{x1D6FF}}$ . It follows that $\unicode[STIX]{x2202}\tilde{\unicode[STIX]{x1D6E4}}/\unicode[STIX]{x2202}\tilde{t}\sim (\tilde{t}/\mathit{Re})^{1/2}\tilde{\unicode[STIX]{x1D6E4}}$ , and scaling time by $\mathit{Re}^{-1/3}$ eliminates dependence on $\mathit{Re}$ . This result also indicates that at early times $\log \tilde{\unicode[STIX]{x1D6E4}}\sim \tilde{t}^{3/2}/\mathit{Re}^{1/2}$ , consistent with a series expansion of equation (3.23) about $\tilde{t}=0$ , which yields

(3.26) $$\begin{eqnarray}\unicode[STIX]{x1D703}(\tilde{t})=1-\frac{4}{3\sqrt{\unicode[STIX]{x03C0}}}\left(\frac{(2\unicode[STIX]{x03C0})^{4/3}}{\mathit{Re}^{1/3}}\tilde{t}\right)^{3/2}\end{eqnarray}$$

to leading order in $\tilde{t}$ .

3.2 Stokes flow limit

In the Stokes flow limit, i.e.  $\mathit{Re}\rightarrow 0$ , of (2.10), the subsurface flow decays exponentially without oscillation (see also Masoud & Shelley Reference Masoud and Shelley2014). In this limit, the Laplace-domain solution for the surfactant evolution simplifies to

(3.27) $$\begin{eqnarray}\breve{\unicode[STIX]{x1D6E4}}(k,s)=\frac{1}{s+\unicode[STIX]{x03C0}}\hat{\unicode[STIX]{x1D6E4}}_{0}(k),\end{eqnarray}$$

and inversion provides

(3.28) $$\begin{eqnarray}\hat{\unicode[STIX]{x1D6E4}}(k,\tilde{t})=\exp (-\unicode[STIX]{x03C0}\tilde{t})\hat{\unicode[STIX]{x1D6E4}}_{0}(k).\end{eqnarray}$$

Furthermore, we obtain

(3.29a ) $$\begin{eqnarray}\displaystyle & \displaystyle \breve{p}(k,{\tilde{y}},s)=\breve{p}_{0}(k,s)\,\exp (-|k|{\tilde{y}}), & \displaystyle\end{eqnarray}$$

(3.29b ) $$\begin{eqnarray}\displaystyle & \displaystyle \breve{u}(k,{\tilde{y}},s)=-\frac{\text{i}}{2k}\breve{p}_{0}(k,s)\,\exp (-|k|{\tilde{y}})\,(|k|{\tilde{y}}-1), & \displaystyle\end{eqnarray}$$

(3.29c ) $$\begin{eqnarray}\displaystyle & \displaystyle \breve{v}(k,{\tilde{y}},s)={\textstyle \frac{1}{2}}\breve{p}_{0}(k,s)\,{\tilde{y}}\exp (-|k|{\tilde{y}}), & \displaystyle\end{eqnarray}$$

$\tilde{\unicode[STIX]{x1D6E4}}_{0}(\tilde{x})=\sin (2\unicode[STIX]{x03C0}\tilde{x})$

(3.30a ) $$\begin{eqnarray}\displaystyle & \displaystyle \tilde{\unicode[STIX]{x1D6E4}}(\tilde{x},\tilde{t})=\sin (2\unicode[STIX]{x03C0}\tilde{x})\exp (-\unicode[STIX]{x03C0}\tilde{t}), & \displaystyle\end{eqnarray}$$

(3.30b ) $$\begin{eqnarray}\displaystyle & \displaystyle \tilde{u} (\tilde{x},{\tilde{y}},\tilde{t})=-{\textstyle \frac{1}{2}}\cos (2\unicode[STIX]{x03C0}\tilde{x})(1-2\unicode[STIX]{x03C0}{\tilde{y}})\exp (-2\unicode[STIX]{x03C0}{\tilde{y}})\exp (-\unicode[STIX]{x03C0}\tilde{t}), & \displaystyle\end{eqnarray}$$

(3.30c ) $$\begin{eqnarray}\displaystyle & \displaystyle \tilde{v}(\tilde{x},{\tilde{y}},\tilde{t})=-{\textstyle \frac{1}{2}}\sin (2\unicode[STIX]{x03C0}\tilde{x})(2\unicode[STIX]{x03C0}{\tilde{y}})\exp (-2\unicode[STIX]{x03C0}{\tilde{y}})\exp (-\unicode[STIX]{x03C0}\tilde{t}), & \displaystyle\end{eqnarray}$$

(3.30d ) $$\begin{eqnarray}\displaystyle & \displaystyle \tilde{p}_{0}(\tilde{x},\tilde{t})=-2\unicode[STIX]{x03C0}\sin (2\unicode[STIX]{x03C0}\tilde{x})\exp (-\unicode[STIX]{x03C0}\tilde{t}). & \displaystyle\end{eqnarray}$$

Figure 4.

Stokes flow solution for the initial surfactant distribution $\tilde{\unicode[STIX]{x1D6E4}}_{0}(\tilde{x})=\sin (2\unicode[STIX]{x03C0}\tilde{x})$ in the subdomain $[0,1]\times [0,1]$ of the semi-infinite domain. The structure of the flow remains constant as it decays exponentially.

3.3 Developing flow under a steady surfactant distribution

To highlight how the time scale of the coupled transient problem differs from the advective and viscous time scales, we consider briefly a simpler problem: the response of the fluid to a sudden and sustained change in surfactant distribution. In contrast, the surfactant distribution in the coupled problem is initially changed and then allowed to evolve over time. Specifically, we consider the response of a motionless fluid to $\tilde{\unicode[STIX]{x1D6E4}}(\tilde{x},\tilde{t})=\sin (2\unicode[STIX]{x03C0}\tilde{x})H(\tilde{t})$ where $H(\cdot )$ is the Heaviside step function. For this condition, we have

(3.31a,b ) $$\begin{eqnarray}\hat{\unicode[STIX]{x1D6E4}}(k,\tilde{t})=-\text{i}\unicode[STIX]{x03C0}\,\text{sgn}(k)\,\unicode[STIX]{x1D6FF}(|k|-2\unicode[STIX]{x03C0})\,H(\tilde{t}),\quad \breve{\unicode[STIX]{x1D6E4}}(k,s)=-\frac{\text{i}\unicode[STIX]{x03C0}}{s}\,\text{sgn}(k)\,\unicode[STIX]{x1D6FF}(|k|-2\unicode[STIX]{x03C0}).\end{eqnarray}$$

Using (3.9a ) and (3.10), the Laplace-domain solution for the surface speed is

(3.32) $$\begin{eqnarray}\breve{u}(k,0,s)=-\frac{4\unicode[STIX]{x03C0}^{3}}{\mathit{Re}\,s^{2}}\left(\sqrt{1+\frac{\mathit{ Re }\,s}{4\unicode[STIX]{x03C0}^{2}}}-1\right)\unicode[STIX]{x1D6FF}(|k|-2\unicode[STIX]{x03C0}),\end{eqnarray}$$

and upon inverting the Laplace (Erdélyi et al. Reference Erdélyi, Magnus, Oberhettinger and Tricomi1954) and Fourier transforms, we obtain the time-domain solution

(3.33) $$\begin{eqnarray}\displaystyle \tilde{u} (\tilde{x},0,\tilde{t})=-\cos (2\unicode[STIX]{x03C0}\tilde{x})\left[\frac{1}{2}\,\text{erf}\left(\sqrt{\tilde{t}^{\prime }}\right)-\tilde{t}^{\prime }\,\text{erfc}\left(\sqrt{\tilde{t}^{\prime }}\right)+\sqrt{\frac{\tilde{t}^{\prime }}{\unicode[STIX]{x03C0}}}\exp (-\tilde{t}^{\prime })\right], & & \displaystyle\end{eqnarray}$$

where $\tilde{t}^{\prime }=4\unicode[STIX]{x03C0}^{2}\tilde{t}/\mathit{Re}$ . This solution is a monotonic rise in time to a steady state. The effect of $\mathit{Re}$ here is a rescaling of time. The time scale for the development of a flow under a steady stress distribution is the viscous time scale $\unicode[STIX]{x1D70C}\unicode[STIX]{x1D706}^{2}/\unicode[STIX]{x1D707}$ in contrast to the advective time scale $t^{\ast }$ of the coupled problem in the Stokes flow limit and $t_{\unicode[STIX]{x1D6E4}}^{\ast }$ in the high $\mathit{Re}$ limit.

4.1 Steady surfactant distribution

As a benchmark, we examine first the transient solution for the flow developing from rest under a steady surfactant distribution, i.e.  $\unicode[STIX]{x1D6E4}(x,t)=\unicode[STIX]{x1D6E4}_{0}(x)~\text{for}~t>0$ . The computed surface speed as a function of time, shown in figure 5 for $\mathit{Re}/(2\unicode[STIX]{x03C0})=0.5$ , 5, and 50 with $\unicode[STIX]{x1D700}=0.1$ , agrees well with the analytical solution, equation (3.33). Though not shown, deviation from the analytical solution begins with $\unicode[STIX]{x1D700}\mathit{Re}/(2\unicode[STIX]{x03C0})\approx 10$ due to the nonlinear term in (2.10) that was neglected in the derivation of the analytical solution.

Figure 5.

Comparison of the analytical solution for the development of the surface flow due to a constant sinusoidal surfactant distribution, equation (3.33), and simulations with $\unicode[STIX]{x1D700}=0.1$ and $\mathit{Re}/(2\unicode[STIX]{x03C0})=0.5$ , 5 and 50. The flow evolves in the same way for all $\tilde{x}$ .

Figure 6.

Numerical and analytical solutions for the decay of small variations in the surfactant concentration with (a) $\mathit{Re}/(2\unicode[STIX]{x03C0})\leqslant 32/27$ and (b) $\mathit{Re}/(2\unicode[STIX]{x03C0})>32/27$ . The dotted curves show the simulation results (with $\unicode[STIX]{x1D700}=0.01$ ), the coloured solid lines show numerical inversion of the Laplace domain solution (3.13), and the solid black lines show the analytical solutions ((3.28) for $\mathit{Re}\rightarrow 0$ and (3.23) for $\mathit{Re}\gg 1$ ). The curves with large dots were computed with sufficiently high $\mathit{Pe}_{s}$ that diffusion is negligible ( $\mathit{Pe}_{s}=1000$ in (a) and in (b) $\mathit{Pe}_{s}=1587$ , 4641 and 10 000 for $\mathit{Re}/(2\unicode[STIX]{x03C0})=4$ , 100 and 1000, respectively), while the curve with smaller dots in (b) shows the computed result for $\mathit{Re}/(2\unicode[STIX]{x03C0})=1000$ with $\mathit{Pe}_{s}=500$ .

4.2 Small $\unicode[STIX]{x1D700}$

Simulations of the full transient problem with decaying small-amplitude ( $\unicode[STIX]{x1D700}\rightarrow 0$ ) variations in the surfactant concentration and high $\mathit{Pe}_{s}$ are shown in figure 6. In the limits of small and large $\mathit{Re}$ , the numerical results approach the corresponding analytical solutions. For intermediate $\mathit{Re}$ , the simulation results were compared with numerical inversion of the Laplace-domain solution (3.14). Since the simulations include the surfactant diffusion term, a sufficiently large value of $\mathit{Pe}_{s}$ is required to keep diffusion negligible. The choice of $\mathit{Pe}_{s}$ must ensure that the diffusion time scale is significantly longer than the surfactant decay time scale, which is not necessarily the advection time scale used to define $\mathit{Pe}_{s}$ . In the limit of low $\mathit{Re}$ , we require

(4.1) $$\begin{eqnarray}\frac{\unicode[STIX]{x1D706}^{2}}{D_{s}}\gg \frac{t^{\ast }}{\unicode[STIX]{x03C0}}\Rightarrow \mathit{Pe}_{s}\gg 1,\end{eqnarray}$$

whereas in the high $\mathit{Re}$ limit the requirement is

(4.2) $$\begin{eqnarray}\frac{\unicode[STIX]{x1D706}^{2}}{D_{s}}\gg \frac{\mathit{Re}^{1/3}t^{\ast }}{(2\unicode[STIX]{x03C0})^{4/3}}\Rightarrow \frac{(2\unicode[STIX]{x03C0})^{4/3}}{\mathit{Re}^{1/3}}\mathit{Pe}_{s}\gg 1,\end{eqnarray}$$

which follows from (3.25). While a fixed high $\mathit{Pe}_{s}$ is suitable at low $\mathit{Re}$ , in which case the surfactant decay time scale is the advection time scale, $\mathit{Pe}_{s}$ must be scaled in proportion to $\mathit{Re}^{1/3}$ for high $\mathit{Re}$ . The simulations in figure 6(a) all had $\mathit{Pe}_{s}=1000$ , while those in figure 6(b) were performed with $\mathit{Pe}_{s}=1587$ , 4641 and 10 000 for $\mathit{Re}/(2\unicode[STIX]{x03C0})=4$ , 100 and 1000, respectively. For these values, $(\mathit{Re}/(2\unicode[STIX]{x03C0}))^{-1/3}\mathit{Pe}_{s}=1000$ . Good agreement in all cases with the Laplace-domain solution for infinite $\mathit{Pe}_{s}$ indicates that these choices of $\mathit{Pe}_{s}$ are satisfactory to suppress the effects of diffusion. For comparison, a simulation with $\mathit{Pe}_{s}=500$ and $\mathit{Re}/(2\unicode[STIX]{x03C0})=1000$ is also shown in figure 6. At this lower $\mathit{Pe}_{s}$ , agreement with the high $\mathit{Pe}_{s}$ solution is still reasonable. Though not shown in figure 6, the surfactant distribution continues to decay faster as $\mathit{Pe}_{s}$ is decreased further. With $\mathit{Pe}_{s}=10$ , the minimum becomes $-0.01$ at $[(2\unicode[STIX]{x03C0})^{4/3}/\mathit{Re}^{1/3}]\tilde{t}\approx 1.5$ .

The evolution of the flow structure at $\mathit{Re}/(2\unicode[STIX]{x03C0})=100$ and $\unicode[STIX]{x1D700}=0.01$ is illustrated in figure 7. As mentioned earlier, the oscillations in the surfactant concentration occur due to the inertia of the flow when the surfactant concentration first becomes uniform. During the initial decay in the amplitude of the surfactant variation, its first zero crossing, and subsequent growth, the flow structure consists of a row of counter-rotating vortices that penetrate deeper into the fluid layer over time. A second layer of vortices forms above the initial layer when the first peak in amplitude (of the opposite sign to the initial condition) is reached. The change in sign of the surface stress therefore reverses the motion of the surfactant and creates this new layer of vortices rather than stopping and reversing the initial vortex layer. The flow structure alternates between having one and two layers of vortices until the third and final peak in amplitude. The flow and surfactant concentration variation then decay monotonically to zero with two rows of vortices present. The development of the flow for $\mathit{Re}/(2\unicode[STIX]{x03C0})=100$ is representative of the structure seen in the oscillatory decay regime ( $\mathit{Re}/(2\unicode[STIX]{x03C0})>32/27$ ). With increasing $\mathit{Re}$ , however, the vortices are compressed vertically.

Figure 7.

Evolution of the flow structure at $\mathit{Re}/(2\unicode[STIX]{x03C0})=100$ and $\unicode[STIX]{x1D700}=0.01$ in the upper half (a), quarter (b–e), and fifth (f) of the simulated domain. The colour indicates the flow speed, and a logarithmic colour scale is used for convenience. Time proceeds from (a) to (f). Different ranges of the vertical axes in (a) and (f) are used to emphasize the structure of the flow.

4.3 Initial conditions with multiple modes

Due to the linearity of the governing equations in the limit of small $\unicode[STIX]{x1D700}$ , the solutions for arbitrary periodic initial surfactant distributions can be constructed by superposition of the solutions for single modes. Expanding $\tilde{\unicode[STIX]{x1D6E4}}_{0}(\tilde{x})$ for $\tilde{x}\in [0,1]$ in a Fourier series as $\sum _{n=1}^{\infty }a_{n}\sin (2\unicode[STIX]{x03C0}n\tilde{x}+\unicode[STIX]{x1D719}_{n})$ , we see that the solution will consist of the single-mode solutions for wavelengths $\unicode[STIX]{x1D706}_{n}=\unicode[STIX]{x1D706}_{1}/n$ where $\unicode[STIX]{x1D706}_{1}$ is the wavelength of the fundamental mode. These modes have corresponding Reynolds numbers $\mathit{Re}_{n}=\mathit{Re}_{1}/n$ and time scales $t_{n}^{\ast }=t_{1}^{\ast }/n$ . The fundamental mode decays slowest, as shown in figure 8 for (a) an initial distribution that has two modes and (b) one that is random. For the simulation with two modes, $\mathit{Pe}_{s}$ was chosen as 20 800 to satisfy equation (4.2) for the $n=3$ mode; for the random distribution, $\mathit{Pe}_{s}=10\,000$ for the fundamental mode. The evolution of the random distribution shows that when the fundamental mode reaches its first zero, the third and higher modes are effectively negligible because the solution at this time is bounded by the amplitude of the second mode.

Figure 8.

The decay of initial surfactant distributions with multiple modes. The time axis is scaled using the Reynolds number for the fundamental mode, $\mathit{Re}_{1}$ , with $\mathit{Re}_{1}/(2\unicode[STIX]{x03C0})=1000$ and $\unicode[STIX]{x1D700}=0.01$ for both cases. In (a), the initial condition is $\tilde{\unicode[STIX]{x1D6E4}}_{0}(\tilde{x})=\cos (2\unicode[STIX]{x03C0}\tilde{x})+(1/4)\sin (6\unicode[STIX]{x03C0}\tilde{x})$ . The numerical solutions (coloured lines with symbols) are compared with numerical inversion of the Laplace-domain solution (solid black lines) at $\tilde{x}=0$ , $5/12$ , $1/2$ and $11/12$ . In (b), $\tilde{\unicode[STIX]{x1D6E4}}_{0}$ is a random value from $[-1,1]$ at each point in the simulated domain. The concentrations are shown at positions $\tilde{x}$ that are multiples of $1/8$ (grey lines with dots), and they are scaled by the magnitude $a_{1}$ of the fundamental mode in the initial condition. The black lines show the decay of the first ( $n=1$ , solid) and second ( $n=2$ , dashed) modes in the initial distribution obtained by numerical inversion of the Laplace-domain solution and evaluated at the positions of maximum amplitude for each mode.

4.4 Finite $\unicode[STIX]{x1D700}$

For small $\unicode[STIX]{x1D700}$ , the portions of the surface with positive and negative deviation from the mean surfactant concentration evolve symmetrically, and this case is described by the analytical and Laplace-domain solutions. The simulations in figure 9 show the asymmetry in the evolution of the surfactant as $\unicode[STIX]{x1D700}$ increases to one, which is its upper bound since negative concentrations are inadmissible. While the results for $\unicode[STIX]{x1D700}=0.01$ match the $\unicode[STIX]{x1D700}\rightarrow 0$ solutions, the region with positive $\tilde{\unicode[STIX]{x1D6E4}}$ initially widens and decays faster for higher $\unicode[STIX]{x1D700}$ with low $\mathit{Re}$ (figure 9 a), the region with negative $\tilde{\unicode[STIX]{x1D6E4}}$ narrows, and the concentration profile sharpens around the minimum. These effects are due to the nonlinear term in (2.9), which accelerates the decrease of $\tilde{\unicode[STIX]{x1D6E4}}$ near its maximum, makes $\unicode[STIX]{x2202}\tilde{\unicode[STIX]{x1D6E4}}/\unicode[STIX]{x2202}\tilde{t}$ non-zero where initially $\tilde{\unicode[STIX]{x1D6E4}}\approx 0$ , and decreases $\unicode[STIX]{x2202}\tilde{\unicode[STIX]{x1D6E4}}/\unicode[STIX]{x2202}\tilde{t}$ near the minimum $\tilde{\unicode[STIX]{x1D6E4}}$ . In fact, as will be described in greater detail for higher $\mathit{Re}$ , the convective term vanishes at the minimum $\tilde{\unicode[STIX]{x1D6E4}}$ as $\unicode[STIX]{x1D700}\rightarrow 1$ leaving only diffusion to drive changes in concentration there.

Figure 9.

Evolution of the surfactant concentration in the limit $\unicode[STIX]{x1D700}\rightarrow 0$ (analytical solutions or numerical inversion of Laplace-domain solution) and numerical solutions for finite-amplitude variation in concentration ( $\unicode[STIX]{x1D700}=0.01$ , 0.5 and 1) at (a) $\mathit{Re}/(2\unicode[STIX]{x03C0})=0.1$ , (b) 100 and (c) 1000.

At the higher $\mathit{Re}$ , for which oscillations occur (figure 9 b,c), the effects of increasing $\unicode[STIX]{x1D700}$ are more complex because of the nonlinear terms in both the surfactant evolution and fluid flow. At early times ( $[(2\unicode[STIX]{x03C0})^{4/3}/\mathit{Re}^{1/3}]\tilde{t}<1$ ), the behaviour is similar to that with lower $\mathit{Re}$ . Unlike at lower $\mathit{Re}$ , the concentration at the minimum then begins to increase faster, and for $\unicode[STIX]{x1D700}=1$ , the growth rate exceeds that in the immediately adjacent regions, and a local maximum and two adjacent minima form in the surfactant concentration. This local maximum is more pronounced and occurs earlier at the highest $\mathit{Re}$ . From $[(2\unicode[STIX]{x03C0})^{4/3}/\mathit{Re}^{1/3}]\tilde{t}\approx 1.5$ to 2, a temporary doubling in the spatial frequency of the concentration profile is evident.

The asymmetry in the evolution of the surfactant concentrations at the locations of the initial minimum and maximum is emphasized in figure 10. The most significant change in the dynamics with increasing $\unicode[STIX]{x1D700}$ is the emergence of a delay in the growth of the concentration from $\tilde{\unicode[STIX]{x1D6E4}}=-1$ . In figure 9, this was seen as a narrowing of the region with negative $\tilde{\unicode[STIX]{x1D6E4}}$ . For $\unicode[STIX]{x1D700}$ approaching sufficiently close to one, the advective terms in the surfactant dynamics become negligible for any $\mathit{Pe}_{s}$ , leaving only diffusion to drive changes in $\tilde{\unicode[STIX]{x1D6E4}}$ at early times. The diffusive flux that determines the duration of the initial delay depends on the initial concentration profile and the developing flow, which sharpens the profile, increasing the flux. The dominance of the diffusive term persists until $\tilde{\unicode[STIX]{x1D6E4}}$ increases sufficiently that the advective terms take over. Due to the effect of diffusion on the initial delay, the surfactant evolution for high $\unicode[STIX]{x1D700}$ is expected to depend on $\mathit{Pe}_{s}$ .

Figure 10.

Asymmetry in the decay of the surfactant concentration from its maximum and minimum values at $\mathit{Re}/(2\unicode[STIX]{x03C0})=1000$ , $\mathit{Pe}_{s}=10\,000$ , and several $\unicode[STIX]{x1D700}$ . The solid black line shows the numerical inversion of the Laplace-domain solution for the limit of $\unicode[STIX]{x1D700}\rightarrow 0$ , equation (3.14).

4.5 Parameters for physical systems

Where do typical physical systems lie in the dimensionless parameter space and how would they behave? Consider an aqueous surfactant solution with a clean surface tension $\unicode[STIX]{x1D6FE}_{0}$ of $72~\text{mN}~\text{m}^{-1}$ and 10 % variation in tension around a mean of $50~\text{mN}~\text{m}^{-1}$ . Using $\unicode[STIX]{x1D70C}=1000~\text{kg}~\text{m}^{-3}$ and $\unicode[STIX]{x1D707}=1~\text{mPa}~\text{s}$ , we obtain $\mathit{Re}=2.2\times 10^{4}$ for a millimetre-scale variation in surface tension. The advective time scale is $45~\unicode[STIX]{x03BC}\text{s}$ . Assuming a diffusivity of the order of $10^{-9}~\text{m}^{2}~\text{s}^{-1}$ , which is typical of small solute molecules, the diffusive time scale is $\unicode[STIX]{x1D706}^{2}/D_{s}=10^{3}~\text{s}$ and $\mathit{Pe}_{s}=2.2\times 10^{7}$ . This system is clearly in the high $\mathit{Re}$ limit and the surfactant decay time scale is $1.3~\text{ms}$ , much shorter than the diffusive time scale. For centimetre-scale variation, $\mathit{Re}=2.2\times 10^{5}$ , the advective time scale is $0.45~\text{ms}$ , the diffusive time scale is $10^{5}~\text{s}$ , $\mathit{Pe}_{s}=2.2\times 10^{8}$ and the surfactant decay time scale is $27~\text{ms}$ . For both these sample sets of parameters, the surfactant decay time scales are significantly shorter than the diffusion time scales, satisfying the key assumption in the modelling that $\mathit{Pe}_{s}$ is large and diffusion can be neglected. The diffusivities of larger solute molecules can be up to several orders of magnitude lower, such as $10^{-10}~\text{m}^{2}~\text{s}^{-1}$ for sodium dodecyl sulphate, and values of $\mathit{Pe}_{s}$ for such solutes would be correspondingly higher. Considering that Schmidt numbers $\mathit{Sc}_{s}=\unicode[STIX]{x1D707}/(\unicode[STIX]{x1D70C}D_{s})=\mathit{Pe}_{s}/\mathit{Re}$ for aqueous systems are typically $O(1000)$ , the assumption of high $\mathit{Pe}_{s}$ is reasonable for the oscillatory surfactant decay regime in which $\mathit{Re}$ exceeds one. The effects of diffusion need to be included for lower $\mathit{Re}$ , however, unless the system has a higher $\mathit{Sc}_{s}$ . The high $\mathit{Re}$ of the sample systems, which exceed those considered in this paper for $\unicode[STIX]{x1D700}>0$ , motivate further study of nonlinear phenomena at high $\mathit{Re}$ and finite $\unicode[STIX]{x1D700}$ , in particular the limit $\unicode[STIX]{x1D700}\rightarrow 1$ .