Axisymmetric spreading of surfactant from a point source

TL;DR

This paper theoretically analyzes the steady axisymmetric flow caused by surfactant spreading from a point source on a fluid interface, distinguishing between adsorbed and dissolved phases and deriving universal flow profiles.

Contribution

The study provides a self-similar analytical model for Marangoni-driven flow due to surfactant spreading, highlighting differences between adsorbed and dissolved phases.

Findings

Flow velocity decays as r^{-3/5} for adsorbed phase

Flow velocity decays as r^{-1} for dissolved phase

Universal flow profiles are identified for both cases

Abstract

Guided by computation, we theoretically calculate the steady flow driven by Marangoni stress due to surfactant introduced on a fluid interface at a constant rate. Two separate extreme cases, where the surfactant dynamics is dominated by the adsorbed phase or the dissolved phase are considered. We focus on the case where the size of the surfactant source is much smaller than the size of the fluid domain, and the resulting Marangoni stress overwhelms viscous forces so that the flow is strongest in a boundary layer close to the interface. We derive the resulting flow in a region much larger than the surfactant source but smaller than the domain size is described by approximating it with a self-similar profile. The radially outward component of fluid velocity decays with the radial distance $r$ as $r^{-3/5}$ when the surfactant spreads in an adsorbed phase, and as $r^{-1}$ when it spreads in a dissolved phase. Universal flow profiles that are independent of the system parameters emerge in both the cases. Three hydrodynamic signatures are identified to distinguish between the two cases and verify the applicability of our analysis with successively stringent tests.

Figures

Equations

$$(ru{)}_r+r{w}_z$$

$$(r{u}^2{)}_r+(rwu{)}_z+r{p}_r$$

$$(ruw{)}_r+(r{w}^2{)}_z+r{p}_z$$

$$w=0\text{ and }\mu u_z={\sigma }_r\text{ at }z=0.$$

$$\sigma ={\sigma }_0-{\Gamma }_2{c}_2,$$

$$2\pi ru(r,0)c_2$$

$$\sigma ={\sigma }_0-{\Gamma }_3{c}_3\text{at}z=0,$$

$$u{c}_{3,r}+w{c}_{3,z}$$

$${\int }_{-\infty }^0\left[2\pi ru(r,z)c_3(r,z)-Dr{c}_{3,r}(r,z)\right]\text{d}z=q_3=\text{constant}\text{for all }r.$$

$$u\propto \frac{K_2^{2/5}{\nu }^{1/5}}{r^{3/5}},\delta \propto \frac{{\nu }^{2/5}{r}^{4/5}}{K_2^{1/5}},\text{ and }\frac{c_2}{q_2}\propto \frac{1}{K_2^{2/5}{\nu }^{1/5}{r}^{2/5}}.$$

$$\xi =\frac{z}{\delta (r)},\text{ where }\delta =\frac{{\nu }^{2/5}{r}^{4/5}}{K_2^{1/5}}$$

$$u{u}_r+w{u}_z=\nu u_{zz},p_z=0,(ur{)}_r+r{w}_z=0,$$

$$u(r,z)=\frac{K_2^{2/5}{\nu }^{1/5}}{r^{3/5}}{f}^{\prime }(\xi ),\text{ and }w(r,z)=-\frac{K_2^{1/5}{\nu }^{3/5}}{5r^{4/5}}\left(6f(\xi )-4\xi f^{\prime }(\xi )\right).$$

$$f^{\prime \prime \prime }(\xi )+\frac{3}{5}{f}^{\prime }(\xi {)}^2+\frac{6}{5}f(\xi )f^{\prime \prime }(\xi )=0.$$

$$u_z(r,z=0)=-\frac{{\Gamma }_2{q}_2}{2\pi \mu }{\left(\frac{1}{ru(r,0)}\right)}_r=-K_2{\left(\frac{1}{ru(r,0)}\right)}_r.$$

$$f(0)=0,f^{\prime \prime }(0)f^{\prime }(0)=\frac{2}{5},\text{and}{f}^{\prime }(-\infty )=0.$$

$$f(\xi )=f_{\infty }+ϵ\left(a_1+b_1\xi +c_1{e}^{-6f_{\infty }\xi /5}\right)+O({ϵ}^2),$$

$${\nabla }^2\varphi =\frac{1}{{\rho }^2}{\left({\rho }^2{\varphi }_{\rho }\right)}_{\rho }+\frac{1}{{\rho }^2\sin \theta }{\left(\sin \theta {\varphi }_{\theta }\right)}_{\theta }=0,$$

$$\frac{1}{\rho }{\varphi }_{\theta }(\rho ,\theta \to {\frac{\pi }{2}}^{-})=w(\rho ,\xi \to -\infty )=-1.2f_{\infty }{K}_2^{1/5}{\nu }^{3/5}{\rho }^{-4/5}.$$

$$(1-s^2)f_o^{\prime \prime }(s)-2s{f}_o^{\prime }(s)+n(n+1)f_o=0.$$

$$\begin{array}{rl}& \varphi =C_{\varphi }{r}^{1/5}{P}_{\frac{1}{5}}(\cos \theta ),v=\frac{1}{5}{C}_{\varphi }{r}^{-4/5}{P}_{\frac{1}{5}}(\cos \theta ),q=-C_{\varphi }{r}^{-4/5}{P}_{\frac{1}{5}}^{\prime }(\cos \theta )\sin \theta ,\\ & \text{where }{C}_{\varphi }=\frac{6f_{\infty }{K}_2^{1/5}{\nu }^{3/5}}{5P_{\frac{1}{5}}^{\prime }(0)}.\end{array}$$

$$p=-\frac{C_{\varphi }^2}{2{\rho }^{8/5}}\left(\frac{P_{\frac{1}{5}}^2(\cos \theta )}{25}+P_{\frac{1}{5}}^{\prime 2}(\cos \theta ){\sin }^2\theta \right).$$

$$\frac{{\nu }^2}{K_2}\ll r.$$

$$r\gg \sqrt{\frac{D_2^5}{K_2^2\nu }}.$$

$$p=-C_p\frac{K_2^{2/5}{\nu }^{6/5}}{r^{8/5}}\text{where}{C}_p=\frac{18f_{\infty }^2}{25}\left(\frac{P_{\frac{1}{5}}^2(0)}{25P_{\frac{1}{5}}^{\prime 2}(0)}+1\right)\approx \mathrm{0.955092285086512.}$$

$$\frac{1}{r}{\left(r{H}_r\right)}_r-\frac{H}{l_c^2}=\frac{p}{\sigma }=-C_p\frac{K_2^{2/5}{\nu }^{6/5}}{\sigma r^{8/5}}.$$

$$H(r)=-C_p\frac{K_2^{2/5}{\nu }^{6/5}{l}_c^{2/5}}{\sigma }\left[K_0\left(\frac{r}{l_c}\right){\int }_0^{r/l_c}\frac{I_0(s)}{s^{3/5}}ds-I_0\left(\frac{r}{l_c}\right){\int }_{\infty }^{r/l_c}\frac{K_0(s)}{s^{3/5}}ds\right],$$

$$H\approx C_p\frac{K_2^{2/5}{\nu }^{6/5}{l}_c^2}{\sigma r^{8/5}}.$$

$$\frac{1}{r}{\left(r{H}_r\right)}_r\approx -C_p\frac{K_2^{2/5}{\nu }^{6/5}}{\sigma r^{8/5}}⟶H\approx C_p\frac{K_2^{2/5}{\nu }^{6/5}}{\sigma }\left(H_0{l}_c^{2/5}-\frac{25r^{2/5}}{4}\right),$$

$$|H_r(r)|\approx C_p\frac{K_2^{2/5}{\nu }^{6/5}}{\sigma }\frac{5}{2r^{3/5}}\ll 1⟶r\gg {\left(\frac{5C_p}{2}\right)}^{5/3}\frac{K_2^{2/3}{\nu }^2}{{\sigma }^{5/3}}.$$

$$v=\nu \frac{\hat{f}(\theta )}{\rho },q=\nu \frac{\hat{g}(\theta )}{\rho },\text{ and }{c}_3=\frac{\nu }{{\Gamma }_3}\frac{\hat{h}(\theta )}{\rho },$$

Full text

1 Introduction

Surfactant spreading on a liquid has received much attention for thin films (Craster & Matar, 2009), but not for deep layers of fluids. Past attempts at analyzing Marangoni-stress driven flow on a deep layer of fluid have been in the context of thermo-capillary (Bratukhin & Maurin, 1967, Napolitano, 1979, Zebib, Homsy & Meiburg, 1985, Carpenter & Homsy, 1990, e.g.) or thermo-soluto-capillary convection (Bratukhin & Maurin, 1968), with the notable exception of Jensen (1995) who analyzed transient dynamics from localized release of adsorbed surfactant. Interest has recently increased in the study of steady flow set by release of soluble amphiphiles at a constant rate (Roché, Li, Griffiths, Le Roux, Cantat, Saint-Jalmes & Stone, 2014, Le Roux, Roché, Cantat & Saint-Jalmes, 2016). In this case, the surfactant is removed from the vicinity of the interface as it dissolves in fluid bath, establishing a state that changes very slowly with the surfactant concentration in the bath. Consistent theoretical treatment of this type of Marangoni-driven surfactant advection is of fundamental interest and the topic of this article.

The general surfactant transport process coupled with the sorption kinetics and driven by a self-imposed Marangoni stress (Dukhin, Kretzschmar & Miller, 1995, Noskov, 1996, Eastoe & Dalton, 2000, Young, Booty, Siegel & Li, 2009, Xu, Booty & Siegel, 2013) is a problem with formidable complexity. The surfactant concentration is governed by the equilibrium isotherm and the transient dynamics of adsorption and desorption. The transport of surfactant is governed by the surfactant diffusion and advection by the flow established by the Marangoni stress. The Marangoni stress itself depends on the relation between the instantaneous surfactant concentration at the interface and the surface tension. The flow that develops is governed by the Navier-Stokes equations, depending on the density and viscosity of the fluid, and any geometric parameters describing the fluid domain. An analytical solution to the general coupled problem, while greatly desired, is not available.

For the case under consideration, a steady release of soluble amphiphilic surfactant through a point source on the interface has been reported to establish a quasi-steady flow near the surface of a deep pool (Roché et al., 2014, Le Roux et al., 2016). The radial extent of this flow is finite for the surfactants studied (Le Roux et al., 2016) and the axi-symmetric radial velocity profile on the surface of the fluid appears to be universal in shape (Roché et al., 2014). In these studies, it was tacitly assumed that the surfactant spreads in a dissolved phase, and that the surface tension gradient driving the flow arises from the bulk concentration of the surfactant near the interface. More recent experiments by Mandre, Akella, Singh, Singh & Bandi (2017) revealed a power-law decay of the surface velocity either with the distance $r$ from the surfactant source as $r^{-3/5}$ or as $r^{-1}$, and an accompanying depth-wise self-similar boundary layer profile. Here we theoretically expound a possible reason for the two different exponents and the self-similar velocity profile in terms of two extremes in the surfactant dynamics. In this article, we show that when the surfactant dynamics are dominated by the surface adsorbed phase, the fluid velocity decays as $r^{-3/5}$, and when dominated by the dissolved phase the decay is as $r^{-1}$. The former corresponds to the surfactant dynamics in the Marangoni regime where the hydrodynamic time scale is so much faster than the sorption kinetics that little exchange of surfactant occurs between the surface and the bulk. The latter corresponds to the Gibbs regime with the sorption kinetics occuring much faster than the hydrodynamic time scale so that an equilibrium between surface and bulk surfactant may be assumed to have established.

Using numerical computations as initial guides, we develop here the mathematical solutions describing the axisymmetric flow resulting from a concentrated steady source of surfactant on the surface of a deep fluid layer. We maintain the two-way coupling between the surfactant transport and the flow. Instead of the more general problem of soluble surfactant dynamics, we consider the two extreme possibilities where the surfactant dynamics are dominated by the adsorbed phase (i.e. Marangoni regime) or the dissolved phase (i.e. Gibb’s regime). The specific mathematical models for the two cases of adsorbtion-dominated or dissolution-dominated surfactant transport are described in §2.1-2.2. Each case is simulated numerically (see §2.3) and the resulting flow in the region much larger than the surfactant source but much smaller than the flow domain is rationalized using similarity solutions (§3.1 and 3.3). In the case of adsorption-dominated surfactant transport, we exploit a thin boundary layer structure of the flow near the fluid surface to make analytical progress. The criteria for validity of the assumptions that underlie the boundary layer similarity solution for this case are also presented in §3.2. The similarity solution for dissolution-dominated surfactant transport is available due to Bratukhin & Maurin (1967), which we specialize to the limit where the flow occurs in a boundary layer. We conclude in §4 by presenting invariant hydrodynamic signatures that distinguish between the two extremes in surfactant dynamics.

2 Mathematical model

Consider an semi-infinite bath of fluid with a free interface along the $z=0$ plane as schematically shown in Figure 1. A surfactant is released at a constant rate from a point source located at the origin. The axisymmetric forcing suggests description of the flow in cylindrical polar coordinates $(r,z)$ using the radial and axial components of velocity $u(r,z)$ and $w(r,z)$ respectively.

The fluid flow satisfies the incompressible Navier-Stokes equations

[TABLE]

where $p$ is the fluid pressure field divided by fluid density. The fluid is quiescent far from the interface. The interface is assumed to be fixed at $z=0$, and the Marangoni stress arising from the non-uniform surface tension $\sigma$ implies

[TABLE]

2.1 Case of adsorbtion-dominated surfactant dynamics

In case the sorption kinetics are not sufficiently fast relative to the hydrodynamic time scale, the surface tension of the interface is dominated by the dynamics of the adsorbed surfactant. The surface tension $\sigma$ depends on the area concentration of the surfactant $c_{2}$, which we approximate for small concentrations to be linear as

[TABLE]

where $\sigma_{0}$ is the interfacial tension without surfactant and $\Gamma_{2}$ is a material constant. We assume the diffusion to be weak compared to advection, as is the case for most surfactants, so that diffusion of the surfactant can be neglected. A quantitative criteria for the validity of this neglect is developed later in the article. The surfactant is transported along the interface by advection implying

[TABLE]

where $q_{2}$ is the strength of the point source. Note that due to the linear relation between $c_{2}$ and $\sigma$ in (3) and between $c_{2}$ and $q_{2}$ in (4), the parameters $q_{2}$ and $\Gamma_{2}$ only influence the flow through the combination $K_{2}=\Gamma_{2}q_{2}/2\pi\mu$. The only independent parameters in the problem are $K_{2}$ and $\nu$, both possessing purely kinematic dimensions. The focus of this paper is on the steady flow that is established far behind the surfactant spreading front or as the surfactant dissolves in the bulk far away from the source and is depleted from the interface. This condition is implemented computationally by introducing a surfactant sink on the outer boundary of the computational domain.

2.2 Case of dissolution-dominated surfactant dynamics

In the case of a soluble surfactant, the volumetric concentration $c_{3}$ dominates surfactant transport and governs the surface tension profile, which we represent as

[TABLE]

where $\Gamma_{3}$ is a material parameter. The surfactant, in this case, is transported within bulk of the fluid by an advection-diffusion process as

[TABLE]

where $D$ is the diffusivity of the surfactant. Note that the bulk diffusivity of the surfactant may not be neglected in this case, no matter how small its may be. It is so because diffusion across the depth of the fluid layer governs the surface concentration, which in turn determines the Marangoni force. Conservation of the surfactant is expressed in terms of the integrated flux of surfactant $q_{3}$ crossing a cylinder of radius $r$ as

[TABLE]

Similar to the case of insoluble surfactant, the linear relations between $\sigma$, $c_{3}$ and $q_{3}$ imply that $\Gamma_{3}$ and $q_{3}$ only appear in the combination $K_{3}=\Gamma_{3}q_{3}/2\pi\mu$. Thus $K_{3}$, $\nu$ and $D$ are the independent parameters describing the problem, all of them possessing purely kinematic dimensions. The Schmidt number is defined as $\text{Sc}=\nu/D$.

Although we have reduced the number of parameters in each case, they are still too numerous to furnish a useful non-dimensionalization. For example, a length scale $\nu^{2}/K_{2}$ and a velocity scale $K_{2}/\nu$ may be constructed for the insoluble case, but as we find later, these parameter combinations do not represent the scales of length and velocity realized in the solution of the governing equations. It is so because, by supposition, the region of interest for our analysis spans distances from the source much greater than $\nu^{2}/K_{2}$. It is not a priori obvious how the length scale $\nu^{2}/K_{2}$ and the distance from the source must be combined to derive the appropriate length scale for the Marangoni flow (see Napoitano, 1979, for a more detailed dimensional analysis). Therefore, we first attempt a computational solution of the dimensional governing equations, and upon examining their structure construct a suitable non-dimensionalization.

2.3 Numerical solutions

Transient versions of the governing equations were implemented and solved using COMSOL’s Computational Fluid Dynamics and Mathematical Modeling module starting from a static fluid layer with a clean interface until a steady state is reached. The domain was chosen to be a cylinder large enough (radius 2.5 m and depth 2.5 m for the insoluble case, and radius 0.48 m and depth 0.48 m for the soluble case) to approximate an infinite domain, but small enough so that for the parameters chosen no flow instabilities appeared. The $z<0$ region of the cylinder is filled with a fluid possessing density in the range $10^{3}-10^{4}$ kg/m3 and dynamic viscosity $\mu$ ranging from $10^{-1}-10^{-3}$ Pa s. The point source of surfactant located at the origin is assigned a small finite extent (radial extent $r_{0}\approx 10^{-3}$ m) so that the singularity at $r=0$ is regularized. Outflow (i.e. zero gage pressure) boundary condition are imposed at the outer radius and lower boundary to better approximate an infinite domain. Imposing other boundary conditions such as no-slip or slip without penetration does not significantly affect the flow in the region much larger than $r_{0}$ but smaller than the domain size. To achieve a steady state, the surfactant is absorbed on the outer boundary of the cylinder by imposition of $c_{2,3}=0$, while a no-flux condition applies everywhere else on the boundary. We used an non-uniform unstructured triangular mesh and a non-uniform structured rectangular mesh for discretization. The discretization is finer near the surfactant source and the interface to resolve the presence of the boundary layer. The grid was successively refined to test for numerical convergence and confirm the self-similar flow structure.

Sample steady state radial flows for four cases are shown in Figure 2 and 3. The flow is fastest at the interface, and falls off rapidly with depth over a length of O(1 cm) indicating a boundary layer structure driven by Marangoni stresses. The dependence of the velocity and length scales on the problem parameters will be explained using a boundary layer approximation and self-similar solution of the governing equations.

A simple description of the flow in terms of power law decay away from the surfactant source may be derived using dominant balances as follows. At a radius $r$ from the source, let us assume that the radial velocity decays as a power law denoted by $u(r)$ and the surfactant concentration as $c_{2,3}(r)$. A balance between inertia (per unit mass) which scales as $u^{2}/r$ and viscous forces which scale as $\nu u/\delta^{2}$, where $\delta(r)$ is the expected boundary layer thickness, implies $\delta=(\nu r/u)^{1/2}$. When the surfactant spreads in an adsorbed phase, its conservation implies $2\pi ruc_{2}=q_{2}$. Combining this with a balance between the Marangoni stress $\Gamma_{2}c_{2}/r$ and the shear stress on the interface $\mu u/\delta$, yields

[TABLE]

A similar scaling analysis for the case of dissolved surfactant may be used to rationalize the $r^{-1}$ scaling of velocity but not the dimensional pre-factor. It can readily be seen from momentum conservation that if $u\propto r^{n}$, then $\delta\propto r^{(1-n)/2}$. Surfactant transport occurs in a layer of thickness $\delta\text{Sc}^{-1/2}\propto r^{(1-n)/2}$. The Marangoni stress at the surface scales as $u/\delta\propto r^{(3n-1)/2}$, which implies the surfactant concentration scales as $c_{3}\propto r^{(3n+1)/2}$. Therefore, the surfactant flux is $q_{3}\propto ruc_{3}\delta\propto r^{2(1+n)}$. Since $q_{3}$ is independent of $r$, we conclude $n=-1$, implying $u\propto r^{-1}$, $\delta\propto r$, and $c_{3}\propto r^{-1}$. The reason for the failure of this simple scaling argument to yield the dimensional pre-factors will be clarified later.

In Figure 4, we plot the numerically obtained radial velocity on the surface as a function of $r$ and compare it with the afore-derived power laws. (The kink in the profiles at $r=10^{-3}$ m for the dissolution-dominated case corresponds to the abrupt jump in the surfactant flux from a non-zero value inside the source to zero outside.) The power-law behavior compares well with the numerical solution in an intermediate range of radial distances, which are much larger than the source but smaller than the domain size.

3 Self-similar profiles

The power-law description paves way for a self-similar description of the flow in this intermediate region, which we derive next. The simple balances presented above, including the dimensional pre-factors for the adsorbed surfactant case, are closely associated with a self-similar flow with a universal velocity and surfactant profile, which is the topic of the rest of the article.

3.1 Adsorption-dominated surfactant transport

The power law in the surface velocity profile implies scale free dynamics, and therefore a self-similar flow. The self-similar flow can be described using a similarity variable

[TABLE]

is the boundary layer thickness according to (8). Assuming $\delta\ll r$, the governing fluid equations simplify to the Prandtl boundary layer equations

[TABLE]

in a region of length $O(\delta)$ near the interface. Outside this region, the flow is weak, and therefore $u$, $w$ and $p$ are neglected there. A separate treatment of the outer region follows later.

Continuity (1c) may be satisfied by using a self similar form for the velocity components with the dimensional pre-factor determined from (8) as

[TABLE]

where $f$ is function to be determined.

Substituting this ansatz in (10) yields

[TABLE]

The boundary conditions (2), (3) and (4) can be combined into one by eliminating $c_{2}(r)$,to get

[TABLE]

Upon substitution of the ansatz (11) in (13), the condition (2) and the stagnation of the fluid outside the boundary layer yields a full set of boundary conditions for $f$ as

[TABLE]

The simple dominant balance, which we used to derive the power-law expressions in (8), can be seen throughout this derivation as the balance between coefficients of the terms representing the respective physical effects. In this manner, the simple dominant balance analysis represents the more detailed derivation based on self-similarity.

Equations (12) and (14) are numerically solved using a shooting method, which we outline next. The method starts with a guess for $f^{\prime\prime}(0)$ and sets the corresponding $f^{\prime}(0)$ using (14). The initial value problem with the guessed initial condition is then solved numerically using a fourth order Runge Kutta method and its asymptotic behavior as $\xi\to-\infty$ is examined. The solution asymptotically either diverges to $\infty$ or $-\infty$ as $f\sim C\xi^{2/3}$, thereby violating the far-field boundary condition in (14). However, between the cases that diverge to $\infty$ and those that diverge to $-\infty$ is one solution that remains bounded. For this case $f$ approaches a constant value, $f_{\infty}$. This possibility may be examined by making the ansatz $f=f_{\infty}+\epsilon g(\xi)+\dots$ for $\epsilon\ll 1$, resulting in

[TABLE]

where $a_{1}$, $b_{1}$ and $c_{1}$ represent arbitrary constants of integration. Note that $b_{1}$ must be zero so that asymptotic ordering of the solution is maintained as $\xi\to-\infty$ and $a_{1}$ may be absorbed into $f_{\infty}$. The objective of the shooting method is to guess the initial value $f^{\prime\prime}(0)$ such that this bounded solution with $b_{1}=0$ is asymptotically achieved. We find this solution by successive bisection of the interval of $f^{\prime\prime}(0)$ with end points that lead to diverging solutions with opposite signs. This bisection was implemented manually to determine that $f^{\prime\prime}(0)=0.402287361293201$ solves (12-14) numerically to 15 digit accuracy. Correspondingly, $f^{\prime}(0)=0.994314110973191$ and $f_{\infty}=-1.13886447085041$, which leads to $u(r,0)=f^{\prime}(0)K_{2}^{2/5}\nu^{1/5}r^{-3/5}$. The axial velocity just outside the boundary layer is $w(r,\xi\to-\infty)=-1.2f_{\infty}K_{2}^{1/5}\nu^{3/5}r^{-4/5}$. The resulting self-similar profile for $u$ and $w$ based on the solution for $f(\xi)$ and its derivatives is shown in Figure 5.

The similarity solution is compared with the four direct computational solution in Figure 5. The reduction of the variability between the profiles resulting from the differences in $K_{2}$, $\nu$, and $r$, and the excellent comparison with the $f(\xi)$ and its derivatives computed as a solution of (12) and (14), verifies the validity of the similarity solution.

Outside the boundary layer, the flow retains its self-similar structure, but with a different scaling. The viscous stresses may be neglected to leading order, and a potential flow driven by volume flux of fluid entraining the boundary layer may be used to describe the resulting flow. The velocity potential, denoted $\phi(\rho,\theta)$ in terms of spherical polar coordinates $(\rho,\theta)$ defined by $r=\rho\sin\theta$ and $z=-\rho\cos\theta$, with the corresponding velocity components $v=\phi_{\rho}$ and $q=\phi_{\theta}/\rho$ respectively, satisfies

[TABLE]

with the matching condition of the fluid flux into the boundary layer,

[TABLE]

A standard solution of Laplace equation in the form $\phi=C_{\phi}r^{n}f_{o}(\cos\theta)$ may be sought. Substituting this form in the boundary condition yields $n=1/5$, and the Legendre differential equation for $f_{o}$ as

[TABLE]

A solution may be written in terms of the regular Legendre function $P_{\frac{1}{5}}$, which solves this differential equation, as

[TABLE]

(Note that $C_{\phi}$ here has SI units of m9/5/s.)

This solution is compared with the computational solution in Figure 6. To isolate the intermediate range of scales much larger than the source and the boundary layer thickness, but much smaller than the computational domain, this comparison is limited to computational points with 0.02 m $<\rho<$ 0.1 m. The raw values of $v$ and $q$ in this region vary by about two orders of magnitude for a fixed $\theta$, however upon accounting for the dimensional factor $C_{\phi}$ and the spatial scaling factor of $\rho^{-4/5}$ the variation reduces to within 10%.

The momentum balance in the outer region reduces to the Bernoulli equation, and yields the pressure variation as $p=-\dfrac{1}{2}(v^{2}+q^{2})$. Using this relation, the pressure may be estimated to leading order to be

[TABLE]

The pressure variation across the boundary layer is negligible due to its thinness, and therefore the pressure everywhere in the fluid is approximated by (20).

3.2 Asymptotic criteria for the validity of self-similar solution

We now revisit the assumptions made in deriving the similarity solution, and where possible derive the quantitative criteria that for their validity. Perhaps the most striking assumption underlying this analysis is simplification of the surfactant transport processes. Since it is the premise of our approach that the velocity field encodes information about the surfactant transport, we do not examine the validity of the assumptions related to surfactant transport. Instead, in this article we merely pursue the logical conclusion of these assumptions, so that they can be compared against experimental measurements. However, many other assumptions were made in the interest of semi-analytical results, and next we verify their validity.

3.2.1 Boundary layer approximation

The application of the boundary layer approximation implies $\partial_{r}\ll\partial_{z}$ and that the outer (radial) flow is much weaker than the one in the boundary layer. The former implies $r\gg\delta$, or equivalently from (9)

[TABLE]

Furthermore, the neglect of the outer radial velocity, which scales as $C_{\phi}/r^{4/5}$ from (19), compared to the boundary layer radial velocity, which scales as ${K_{2}^{2/5}\nu^{1/5}}/{r^{3/5}}$ from (11), is equivalent to (21). The criterion (21) implies the similarity solution to be valid at radial distances much larger than the length scale that can be constructed from dimensional analysis.

3.2.2 Neglect of surfactant diffusion

According to (4), the surfactant concentration profile also decays as a power law, implying that the radial length scale for variation in $c_{2}$ is $r$. The domination of the advective flux $uc_{2}$, which scales as ${K_{2}^{2/5}\nu^{1/5}}/{r^{3/5}}\times c_{2}$ from (11), over the diffusive flux, which scales as $D_{2}/r\times c_{2}$, $D_{2}$ being the diffusivity of the surfactant on the surface, implies

[TABLE]

3.2.3 Flatness of the interface

We had neglected the interface deformation motivated by the observations of Roché et al. (2014) and Mandre et al. (2017), who report no perceptible deformation of the interface was observed, except for a liquid bridge conecting the interface to the conduit conveying the surfactant solution. The closed form similarity solution allows us to quantitatively estimate the deformation of the interface due to the non-uniform pressure at the interface. Noting that $P_{\frac{1}{5}}(0)=\pi^{1/2}/\Gamma(11/10)\Gamma(2/5)$ and $P^{\prime}_{\frac{1}{5}}(0)=-2\pi^{1/2}/\Gamma(3/5)\Gamma(-1/10)$ (see Abramowitz & Stegun, 1964, §8.6), where $\Gamma$ is the Gamma function, the interface pressure may be approximated as

[TABLE]

To leading order of the flat-interface approximation, a combination of gravity and surface tension maintain the interface flat against this pressure. The change in elevation of the interface, $H(r)$, caused by the fluid pressure is given by

[TABLE]

where $l_{c}$ is the capillary length and $\sigma$ is the surface tension divided by fluid density. While a solution to this equation may formally be written down in terms of Bessel and Lommel functions,

[TABLE]

where the lower limits of integration are chosen to satisfy the conditions of regularity of $H$ at $r=0$ and $r\to\infty$. For $r\gg l_{c}$, the gravitational force dominates over the surface tension, and the interface may be approximated by

[TABLE]

For $r\ll l_{c}$, the Laplace pressure from surface tension dominates over gravity, leading to the dominant balance

[TABLE]

where $H_{0}\approx 6.9484928934459429$ is a constant of integration determined from a numerical evaluation of the integral in (25). The slope of the interface must remain small for the validity of the approximations underlying the similarity solution. This criterion translates to

[TABLE]

3.3 Dissolution-dominated surfactant transport

The boundary layer approximation, that so successfully describes the flow resulting in the case of surfactant transported in the adsorbed phase, fails for the case of dissolved surfactant. The reason for this failure will be presented later. However, a self-similar solution is still possible, as first demonstrated by Bratukhin & Maurin (1967). Here we present a comparison of their solution with our numerical results. In addition, we derive the asymptotic behavior of this solution in the physically relevant limit of $\text{Sc}=\nu/D\gg 1$ and show that a universal self-similar profiles exist in the limit $K_{3}\text{Sc}/\nu^{2}\gg 1$.

The solution is best presented in terms of spherical polar coordinates $(\rho,\theta)$ defined by $r=\rho\sin\theta$ and $z=-\rho\cos\theta$, with the corresponding velocity components $v$ and $q$ repectively. The self-similar ansatz in this case is

[TABLE]

in terms of the symbols $\hat{f}$, $\hat{g}$ and $\hat{h}$. An exact solution of the Navier-Stokes equations is possible in this case, represented in terms of a single parameter $A$ as

[TABLE]

The parameter $A$ quantifies the strength of the Marangoni forcing (the surface velocity is $A\nu/2r$) and is related to the surfactant release rate by

[TABLE]

The results of the numerical solutions shown in Figure 3, when rescaled according to the self-similar ansatz, collapse on the similarity solution (30a), as shown in Figure 7. Note that, while the flow is self-similar, the profile shape is not universal, but varies with $A$, which depends on the strength of Marangoni stress as represented by ${K_{3}\text{Sc}}/{\nu^{2}}$ and Sc. This dependence is shown in Figure 7(d), for Sc from $\sqrt{10}$ to $10^{6}$, the range of Sc that encompasses many common chemical species. This range of Sc corresponds to the asymptotic limit $\text{Sc}\gg 1$. This dependence from (31) has the asymptotic limits

[TABLE]

In these two extremes, we expect a universal flow profile to emerge. In the case we are interested where $A\gg 1$, the solution (30a) simplifies in a region of thickness $\Delta nz=O(1)$ (transforming to cylindrical polar coordinates) to

[TABLE]

where $\hat{\xi}={\Delta nz}/{2r}$, $A$ and $\Delta n=\sqrt{1+A}$ are determined from (31) (or Figure 7(d)). In the limit of large $A$, the profile approaches

[TABLE]

Varying $A$ in (34) simply rescales the boundary layer thickness as $\delta=2r/A^{1/2}$ and the velocity magnitude as $A\nu/2r$, but does not change the leading-order shape of the profile. In this sense, the profile is universal.

The flow profile in (34) in fact describes a free radial jet forced by a steady point momentum source at the origin and a shear free boundary condition on the surface, as originally derived by Squire (1955). Consider the solution of (10), which decays as $1/r$ with the ansatz

[TABLE]

where $a$ is a constant to be determined as part of the solution. The self-similar form $\tilde{f}(\xi)$ satisfies

[TABLE]

The solution to this equation that satisfies $w(r,z=0)=0$ and $u(r,z\to\-\infty)=0$ is

[TABLE]

which is identical to the leading order of (34) with $a=A/2$. This boundary layer approximation applies in the limit $a\gg 1$, equivalent to the one under which (34) was derived. However, note that (37) and (34) are unable to satisfy the Marangoni stress condition at $\xi=0$. Therefore, the value of $a$ remains undetermined at this order of boundary layer theory. Bratukhin and Maurin’s solution is a perturbation of Squire’s jet, where the Marangoni stress is supported by the next term of $O(\Delta n^{-1})$, and therefore the surface shear stress does not scale with dimensional parameters according to the leading order. Due to this structure of the flow, the simple scaling analysis of the type presented in §2.3 fails to capture the dependence of dimensional pre-factors on the parameters. Using the next order correction in (33) or the complete solution in (30a), we see that the shear at the surface is $A\nu/r^{2}$ (weaker than the dimensionally expected $A\nu/\delta r$ as explained). The shear stress on the surface $\mu A\nu/r^{2}$ thus balances the Marangoni stress, which scales as $\Gamma_{3}c_{3}/r$, yielding the scale for surfactant concentration $c_{3}\propto\mu A\nu/\Gamma_{3}r$. Surfactant transport rate $q_{3}$ then scales as $ruc_{3}\delta\text{Sc}^{-1/2}=A^{3/2}\nu^{2}\mu/\Gamma_{3}\text{Sc}^{1/2}$, where we have used that the surfactant boundary layer is thinner than the momentum boundary layer by a factor $\text{Sc}^{1/2}$. This final balance yields $A\propto(K_{3}^{2}\text{Sc}/\nu^{4})^{1/3}$ in agreement with (32) and the results in Figure 7(d).

4 Discussion and conclusion

Here we have presented the similarity solutions underlying the Marangoni stress driven flow due to a steady point source of surfactant, as it spreads either in an adsorbed or dissolved state. The scope of this paper is limited to the case where the flow occurs in a thin layer close to the free surface, as characterized by $\delta\ll r$ for both cases, which in either case results in a universal velocity profile.

Based on the analysis presented in this article, the following hydrodynamic signatures distinguish the two limits, expressed in terms of quantities that may be measured experimentally. The power law exponent for the decay of surface radial velocity $u(r,z=0)$ is the first and obvious signature that distinguishes between the two limits. The surface radial velocity decays proportional to $r^{-3/5}$ in the adsorbtion-dominated case and to $r^{-1}$ in the dissolution-dominated case. The boundary layer scalings that collapse the depth-wise velocity profile constitute the second hydrodynamic signature that distinguish between the two limits. Along with rescaling the velocity $u(r,z)$ with $u(r,z=0)$, the depth-wise coordinates need to be rescaled with factors of $\delta(r)=\sqrt{\nu rf^{\prime}(0)/u(r,z=0)}$ for the profiles to collapse with $f(\xi)/f(0)$ in the adsorption-dominated case. The appropriate depth-wise coordinate rescaling factor can be derived in the dissolution-dominated regime from (34) to be $\sqrt{\nu r/u(r,z=0)}$. The third signature is the relation between the surface velocity and the surface shear stress at different radii. In the adsorption-dominated case, this relation is

[TABLE]

whereas for the dissolution-dominated case it is

[TABLE]

Depending on the available experimental accuracy, successively stricter comparison of the measured velocity profile with the theory presented in this manuscript can be made using these three signatures.

Our results rationalize experimental observations by Mandre et al. (2017) of the power-law decay and boundary layer structure of Marangoni-driven flow. The agreement between the experimental measurements and our theory suggests that the more general surfactant dynamics may under certain circumstances be reduced to simple models akin to Marangoni and Gibbs elasticity of surfactant-laden liquid interfaces. The solutions developed here could also provide insight into the surfactant dynamics based on flow velocimetry in experimental systems, such as by Roché et al. (2014) and Le Roux et al. (2016). The precise criteria under which simplification is possible and conditions for transition are left to be undertaken in the future.