1 Introduction
Self-propelled microdevices powered by converting chemical energy into mechanical work bring to bear functionalities that resemble those of microorganisms. These miniature robots are thus promising candidates for targeted drug delivery, microsurgery, self-assembly and microfluidic manipulation (Sitti Reference Sitti2009; Sitti et al. Reference Sitti, Ceylan, Hu, Giltinan, Turan, Yim and Diller2015). Thanks to the recent advances in micro/nano technology, researchers have designed and tested a wide variety of swimming microrobots, many of them biomimetic, to carry out prescribed tasks in microfluidic systems (Dreyfus et al. Reference Dreyfus, Baudry, Roper, Fermigier, Stone and Bibette2005; Masoud & Alexeev Reference Masoud and Alexeev2010; Masoud, Bingham & Alexeev Reference Masoud, Bingham and Alexeev2012; Wang et al. Reference Wang, Duan, Ahmed, Mallouk and Sen2013; Cheang et al. Reference Cheang, Meshkati, Kim, Kim and Fu2014; Duan et al. Reference Duan, Wang, Das, Yadav, Mallouk and Sen2015; Maggi et al. Reference Maggi, Saglimbeni, Dipalo, De Angelis and Di Leonardo2015; Sánchez, Soler & Katuri Reference Sánchez, Soler and Katuri2015; Dey et al. Reference Dey, Wong, Altemose and Sen2016; Domínguez et al. Reference Domínguez, Malgaretti, Popescu and Dietrich2016). However, despite recent notable progress in creating microswimmers, the development of microscopic machines capable of surfing on fluidic interfaces is still in its infancy. These tiny surfing robots can potentially execute missions that are currently very difficult, or even impossible, to accomplish.
The motion of a microsurfer is typically driven by a gradient of surface tension, which is prompted by a form of activity taking place at the surface of the surfer (Rayleigh Reference Rayleigh1889; Nakata et al. Reference Nakata, Iguchi, Ose, Kuboyama, Ishii and Yoshikawa1997; Lauga & Davis Reference Lauga and Davis2012; Zhao & Pumera Reference Zhao and Pumera2012; Zhang et al. Reference Zhang, Duan, Liu and Sen2013; Masoud & Shelley Reference Masoud and Shelley2014; Masoud & Stone Reference Masoud and Stone2014; Pimienta & Antoine Reference Pimienta and Antoine2014; Würger Reference Würger2014; Girot et al. Reference Girot, Danné, Würger, Bickel, Ren, Loudet and Pouligny2016; Malgaretti, Popescu & Dietrich Reference Malgaretti, Popescu and Dietrich2016). This so-called Marangoni propulsion (Bush & Hu Reference Bush and Hu2006; Lauga & Davis Reference Lauga and Davis2012; Masoud & Stone Reference Masoud and Stone2014) appears in biological systems as well. For instance, there are water-walking insects that secrete surface-active materials to boost their walking speed in emergency situations (Bush & Hu Reference Bush and Hu2006). A major step towards designing advanced Marangoni-driven microrobots is to establish the relation between the prescribed surface activity and the propulsion velocity in various configurations. A configuration of practical interest involves surfing in the presence of a confining solid boundary.
Here, we theoretically examine the surfing motion of particles located at a flat liquid–gas interface atop a layer of fluid, which is bounded by an impermeable solid wall. We consider chemically and thermally active particles that are, respectively, sources of an insoluble chemical species and heat. The self-induced scalar fields (i.e. concentration of species and temperature fields) locally change the surface tension and the consequent gradients in the surface tension lead to the propulsion of the particle. It is natural to surmise that Marangoni surfers are pulled towards the higher surface tension direction. Contrary to this common instinct, we find that the surfers may propel in the lower surface tension direction depending on their geometry and proximity to the solid wall. Specifically, we show that spherical particles propel in the reverse direction under the confinement whereas disk-shaped surfers always move in the anticipated direction. We rationalize these behaviours in terms of the competition between the effects of the asymmetric surface tension acting along the three-phase contact line and the Marangoni flow engendered by the surface tension gradients. In the following, we first derive closed-form expressions for the translational velocity of spheres and disks. We then extend the results by presenting an approximate expression for the propulsion speed of oblate spheroidal particles based on the speeds of spheres and disks. Overall, our findings pave the way for designing microsurfers capable of operating in bounded environments.
2 Problem formulation and solution
Consider a solid particle located at a flat surface, at
$z=0$
, sitting above a layer of Newtonian fluid, with constant viscosity
$\unicode[STIX]{x1D707}$
, that is bounded by an impermeable solid wall at
$z=-H$
(see figure 1). The particle translates with the velocity
$\boldsymbol{U}=U\boldsymbol{e}$
due to a non-uniformity in the surface tension
$\unicode[STIX]{x1D6FE}$
arising either from the discharge of an insoluble agent (e.g. surfactant) or from the release of heat by the particle. Here,
$\boldsymbol{e}$
is the unit vector in the direction of motion.

Figure 1. Marangoni propulsion of (a) spherical and (b) disk-shaped particles at a flat interface above a liquid layer of finite depth. The colour maps at the interfaces represent the surface tension distribution stemming either from the discharge of an insoluble agent or from the release of heat by the particle.
$S_{i}$
and
$l_{p}$
denote, respectively, the liquid–gas interface and three-phase contact line pinned to the particle at a
$90^{\circ }$
contact angle. The translation towards the higher surface tension direction is regarded as the forward propulsion.
To reduce the complexity of the analyses, the following (justifiable) assumptions are made: (i) the liquid–gas interface remains flat over distances much greater than the size of the particle, (ii) the three-phase contact line is pinned at a
$90^{\circ }$
contact angle, (iii) the two-dimensional transport of chemical species and the three-dimensional transport of heat are both dominated by diffusion, (iv) gravitational and inertial effects are negligible, (v) the relevant scalars (i.e. concentration of species and temperature) have linear relations with the surface tension and (vi) the release of the chemical agent and heat are symmetric about the direction of motion, which ensures that the motion of the particle is restricted to translation.
A typical experimental set-up for which the assumptions are valid is an aqueous system with particles of diameter
$D\sim 1\,\unicode[STIX]{x03BC}\text{m}$
and propulsion speed
$U\sim 1\,\unicode[STIX]{x03BC}\text{m}~\text{s}^{-1}$
that produce a chemical species with a diffusion constant
${\mathcal{D}}\sim 10^{-10}\,\text{m}^{2}~\text{s}^{-1}$
. The capillary, Péclet and Reynolds numbers of this system are, respectively,
$Ca=\unicode[STIX]{x1D707}U/\unicode[STIX]{x1D6FE}_{0}\sim 10^{-8}$
,
$Pe=DU/{\mathcal{D}}\sim 10^{-2}$
and
$Re=\unicode[STIX]{x1D70C}DU/\unicode[STIX]{x1D707}\sim 10^{-6}$
. The assumptions are similarly valid for an aqueous system with optically heated particles of diameter
$D\sim 10\,\unicode[STIX]{x03BC}\text{m}$
whose propulsion speed is
$U\sim 10\,\unicode[STIX]{x03BC}\text{m}~\text{s}^{-1}$
. Here, the corresponding dimensionless numbers are
$Ca\sim 10^{-7}$
,
$Pe=DU/\unicode[STIX]{x1D6FC}\sim 10^{-3}$
, and
$Re\sim 10^{-4}$
. The non-dimensional numbers are calculated based on the water properties at room condition (i.e. density
$\unicode[STIX]{x1D70C}=997.0479~\text{kg}~\text{m}^{-3}$
, viscosity
$\unicode[STIX]{x1D707}=8.90\times 10^{-4}\,\text{N}~\text{s}~\text{m}^{-2}$
, thermal diffusivity
$\unicode[STIX]{x1D6FC}=0.143\times 10^{-6}\,\text{m}^{2}~\text{s}^{-1}$
and surface tension
$\unicode[STIX]{x1D6FE}_{0}=71.97\times 10^{-3}\,\text{N}~\text{m}^{-1}$
). In both set-ups, the surface tension is strong enough to keep the interface flat (
$Ca\ll 1$
), the primary mode of heat and mass transfer is diffusion (
$Pe\ll 1$
) and the flow is creeping (
$Re\ll 1$
).
Suppose
$\boldsymbol{u}$
and
$\unicode[STIX]{x1D748}$
are, respectively, the velocity and stress fields, in
$z\leqslant 0$
, corresponding to the Marangoni-driven motion of the particle. Let
$\hat{\boldsymbol{u}}$
and
$\hat{\unicode[STIX]{x1D748}}$
denote, respectively, the velocity and stress fields corresponding to the translation with the velocity
$\hat{\boldsymbol{U}}=\hat{U} \boldsymbol{e}$
of an identical particle at an interface with no surface tension gradients. These fields are free of contact-line singularities, as the particle motion is limited to pure translation (O’Neill, Ranger & Brenner Reference O’Neill, Ranger and Brenner1986; Pozrikidis Reference Pozrikidis2007). Hence, according to the Lorentz reciprocal theorem (see Happel & Brenner Reference Happel and Brenner1983),

where
$S_{p}$
is the wetted area of the particle,
$S_{i}$
represents the liquid–gas interface (see figure 1) and
$\boldsymbol{n}$
is the unit vector outward normal to
$S_{p}$
and
$S_{i}$
(see also Masoud & Stone Reference Masoud and Stone2014). Integrals over bounding surfaces at infinity are zero since velocities decay at least as fast as the inverse distance in the far field. Integrals over the solid walls also vanish due to the zero-velocity boundary condition.
Owing to the no-slip condition,
$\boldsymbol{u}=U\boldsymbol{e}$
and
$\hat{\boldsymbol{u}}=\hat{U} \boldsymbol{e}$
on
$S_{p}$
. Also, the components of
$\boldsymbol{u}$
and
$\hat{\boldsymbol{u}}$
normal to the interface are zero as no fluid exchange takes place there. The balance of shear stress at the interface requires
$(\boldsymbol{n}\boldsymbol{\cdot }\unicode[STIX]{x1D748})\boldsymbol{\cdot }\hat{\boldsymbol{u}}=-\unicode[STIX]{x1D735}_{s}\unicode[STIX]{x1D6FE}\boldsymbol{\cdot }\hat{\boldsymbol{u}}$
and
$(\boldsymbol{n}\boldsymbol{\cdot }\hat{\unicode[STIX]{x1D748}})\boldsymbol{\cdot }\boldsymbol{u}=0$
on
$S_{i}$
, where
$\unicode[STIX]{x1D735}_{s}$
is the surface gradient operator. Since no net external force is applied on the particle, the viscous force
$\int _{S_{p}}\boldsymbol{n}\boldsymbol{\cdot }\unicode[STIX]{x1D748}\,\text{d}S$
exerted on
$S_{p}$
is balanced by the surface tension force
$F_{st}\boldsymbol{e}=\int _{\ell _{p}}\unicode[STIX]{x1D6FE}\boldsymbol{t}\;\text{d}\ell$
acting along the three-phase contact line
$\ell _{p}$
, where
$\boldsymbol{t}$
is the unit vector tangent to
$S_{i}$
and normal to
$\ell _{p}$
. Taking the above relations into account, (2.1) reduces to

where
$-\hat{F}_{d}\boldsymbol{e}=\int _{S_{p}}\boldsymbol{n}\boldsymbol{\cdot }\hat{\unicode[STIX]{x1D748}}\,\text{d}S$
is the fluid drag that the translating particle would have experienced in response to its motion at the interface had there been no Marangoni effect and
$F_{M}=\int _{S_{i}}(\hat{\boldsymbol{u}}/\hat{U} )\boldsymbol{\cdot }\unicode[STIX]{x1D735}_{s}\unicode[STIX]{x1D6FE}\,\text{d}S$
is the contribution of the Marangoni flow.
In the conventional approach, calculating
$U$
entails solving the Stokes flow equations for
$\boldsymbol{u}$
subject to the no-slip and Marangoni stress boundary conditions (e.g. Lauga & Davis Reference Lauga and Davis2012; Würger Reference Würger2014). The use of the reciprocal theorem, instead, allows us to obtain
$U$
by just having
$\hat{\boldsymbol{u}}$
(see (2.2)), which is analytically easier to calculate. Conveniently,
$\hat{\boldsymbol{u}}$
is already derived for the creeping motion of spherical and disk-shaped particles between parallel plates (Ganatos, Pfeffer & Weinbaum Reference Ganatos, Pfeffer and Weinbaum1980; Davis Reference Davis1991). Note that when the particle is located at an equal distance from the walls, due to the symmetry, the velocity fields obtained by Ganatos et al. (Reference Ganatos, Pfeffer and Weinbaum1980) and Davis (Reference Davis1991) are identical to those generated by the translation of the same particle along a flat liquid–gas interface sitting above a liquid layer bounded by a solid boundary. This argument holds when the capillary number is very low, which indicates that the deflection of the interface due to the pressure distribution along the free surface can be ignored. The
$Ca\ll 1$
condition is a stronger constraint when the fluid layer is shallow.
Without loss of generality, hereafter, we set the direction of motion to
$\boldsymbol{e}=\boldsymbol{e}_{x}$
(see figure 1). Hence, following the solution of Ganatos et al. (Reference Ganatos, Pfeffer and Weinbaum1980) for spheres of radius
$R$
, the in-plane components of velocity
$\hat{\boldsymbol{u}}$
in the spherical coordinates
$(r,\unicode[STIX]{x1D703},\unicode[STIX]{x1D719})$
at the interface (
$\unicode[STIX]{x1D703}=0$
) are




















Following our assumption, the surface tension changes linearly with the concentration of the surface-active agent and the interface temperature as

where
$\unicode[STIX]{x1D6FE}_{0}$
,
$K_{C}$
and
$K_{T}$
are constants (Acree Reference Acree1984; Adamson Reference Adamson1990). Here,
$C$
and
$T$
denote the quasi-steady average-subtracted concentration and temperature fields that satisfy
$\unicode[STIX]{x1D6FB}_{s}^{2}C=0$
and
$\unicode[STIX]{x1D6FB}^{2}T=0$
, respectively. This is true when the dominant mechanism of transport is diffusion (i.e.
$Pe\ll 1$
). Otherwise, advection–diffusion equations for the concentration and temperature fields must be solved simultaneously along with the Stokes equations for the velocity field
$\boldsymbol{u}$
, in which case the reciprocal theorem no longer provides a shortcut solution.
Inspection of (2.3) and (2.4) reveals that only the first harmonic modes of
$\unicode[STIX]{x1D719}$
appear in
$\hat{\boldsymbol{u}}$
. This, in conjunction with (2.5), indicates that just the dipolar terms of
$C$
and
$T$
contribute to the integral in (2.2). Retaining the dipolar parts of their multipole expansions, we can write the concentration and temperature distributions at the interface as






3 Exact results for spheres and disks
Substituting (2.3)–(2.6) into (2.2), we obtain the propulsion speed for chemically and thermally active spheres and disks. We normalize
$U$
by its corresponding value for a semi-infinite fluid layer
$U_{\infty }$
and plot the results as a function of
$H/R$
(see figure 2). When the fluid layer is very deep (
$H/R\geqslant 100$
) the speed is within
$2\,\%$
of its asymptotic value for all four cases considered, among which chemically and thermally active spheres have, respectively, the slowest and fastest approach rates to
$U_{\infty }$
. On the other hand, when the fluid layer is shallow (
$H/R<1.5$
), we see distinctly different behaviours for spheres and disks. Notably, we find that spheres come into a halt at a certain distance from the wall (
$H/R\simeq 1.25$
), closer than which they propel in the lower surface tension direction (see figure 2
a). This contradicts the default notion that Marangoni surfers self-propel towards the direction with a higher surface tension.

Figure 2. The normalized propulsion speed
$U/U_{\infty }$
as a function of the thickness of the fluid layer
$H/R$
for (a) spherical and (b) disk-shaped particles. The squares and circles represent the results for chemically and thermally active particles, respectively. The insets in (a) and (b) show the variation of
$\hat{F}_{d}/\hat{F}_{d,\infty }$
(normalized Stokes drag experienced by the particle translating at an interface with no surface tension gradients) versus
$H/R$
obtained from the solutions of Ganatos et al. (Reference Ganatos, Pfeffer and Weinbaum1980) and Davis (Reference Davis1991). Here, the speed and drag corresponding to a semi-infinite fluid layer (
$H/R\rightarrow \infty$
) are denoted, respectively, by
$U_{\infty }$
and
$\hat{F}_{d,\infty }$
, whose values for spheres and disks are shown in figure 3. Also, the thin solid line in (a) corresponds to
$U=0$
.
According to (2.2), the sign of
$U$
is set by the competition between two opposing influences, namely the net surface tension force acting along the three-phase contact line and the Marangoni flow induced by the surface tension gradients that are represented by
$F_{st}$
and
$F_{M}$
, respectively. Note that
$\hat{F}_{d}/\hat{U}$
is always positive (see insets in figure 2). Thus, if the contribution of
$F_{M}$
is neglected,
$U$
always has the same sign as
$F_{st}$
, meaning that the particles are always pulled towards the higher surface tension direction. This could be the origin of the common perception about the direction of the Marangoni propulsion.
The effect of the confinement on the propulsion speed of disks is also noteworthy. Figure 2(b) shows that
$U$
for chemically active disks approaches zero as
$H/R\rightarrow 0$
, whereas, in the same limit, the speed of thermally active disks increases logarithmically. The decay of
$U$
in the former can be attributed to the linear rise of
$\hat{F}_{d}$
due to the confinement (see inset in figure 2
b). The speed surge in the latter however stems from the fact that the solid wall is insulated. The combination of adiabatic condition and narrow layer results in steeper temperature gradients, which further augment
$F_{st}$
over
$F_{M}$
. Ultimately, the overall gain of (
$F_{st}+F_{M}$
) outweighs the increase of
$\hat{F}_{d}$
(see figure 2
b and (2.2)).
Thanks to the reciprocal theorem, we have so far analytically calculated
$U$
for spheres and disks. Obtaining
$U$
in a similar manner for general oblate spheroids requires a closed-form expression for
$\hat{\boldsymbol{u}}$
, which is not readily available. In what follows, we propose an approximate method for calculating the propulsion speed of oblate spheroids based on our results for spheres and disks.
4 Approximate results for oblate spheroids
Consider an oblate spheroid translating at an interface above a semi-infinite layer of fluid (
$H/R\rightarrow \infty$
). If we plot the analytically known drag
$\hat{F}_{d,\infty }$
as a function of the particle aspect ratio
$\unicode[STIX]{x1D700}$
(see figure 3
a), we observe that the drag changes almost linearly with the aspect ratio (see also the Appendix of Stone & Masoud (Reference Stone and Masoud2015)). In fact, the curve can be approximated extremely well by

where the superscripts
$S$
and
$D$
denote the values for sphere (
$\unicode[STIX]{x1D700}=1$
) and disk (
$\unicode[STIX]{x1D700}=0$
), respectively. The inset in figure 3(a) shows that the maximum relative error is
$1.3\,\%$
belonging to
$\unicode[STIX]{x1D700}=0.4$
. That the particle resistance varies quasi-linearly with the aspect ratio suggests that the particle mobility (e.g.
$U$
) changes inverse linearly with
$\unicode[STIX]{x1D700}$
. Remarkably, we find that the curves of
$U_{\infty }$
versus
$\unicode[STIX]{x1D700}$
for chemically and thermally active spheroids are closely represented by (see figures 3
b and 3
c)

Again, the maximum error occurs at
$\unicode[STIX]{x1D700}=0.4$
and its magnitude is
$1.1\,\%$
and
$12.5\,\%$
for chemically and thermally active particles, respectively (see insets of figures 3
b and 3
c).

Figure 3. Comparison between the exact and approximate results for oblate spheroidal particles translating at an interface above a semi-infinite layer of fluid (
$H/R\rightarrow \infty$
). (a) The variation of non-dimensional drag
$\hat{F}_{d,\infty }/\unicode[STIX]{x1D707}\hat{U} R$
versus the particle aspect ratio
$\unicode[STIX]{x1D700}$
. The dimensionless propulsion speed of (b) chemically and (c) thermally active spheroids (respectively,
$U_{\infty }\unicode[STIX]{x1D707}/K_{C}B_{C}$
and
$U_{\infty }\unicode[STIX]{x1D707}/K_{T}B_{T}$
) as a function of
$\unicode[STIX]{x1D700}$
. The solid lines represent exact analytical results whereas square symbols depict approximate values obtained via interpolations between the results for spheres and disks. The insets show the relative error of the approximations.
These observations strongly recommend employing the same type of interpolations used in (4.1) and (4.2) for approximating the drag
$\hat{F}_{d}$
and propulsion speed
$U$
of oblate spheroidal particles in the presence of the confining solid boundary. To test the validity of such estimations, we numerically calculate
$\hat{F}_{d}$
and
$U$
for
$\unicode[STIX]{x1D700}=0.4$
and compare the results with those obtained by interpolating between the data presented in figure 2. Our numerical approach integrates a boundary singularity method with the (infinite) image technique (Blake Reference Blake1971). Also, we choose
$\unicode[STIX]{x1D700}=0.4$
because, very likely, the error is maximum at this aspect ratio (see insets of figure 3). Indeed, we find that the approximations for
$U$
are reasonably accurate, particularly for
$H/R\geqslant 4$
(see figure 4). The interpolation performs even better for
$\hat{F}_{d}$
, where the error is below
$5\,\%$
for
$H/R\geqslant 2$
.

Figure 4. Comparison between the numerical and approximate results for oblate spheroidal particles of aspect ratio
$\unicode[STIX]{x1D700}=0.4$
. The squares and circles, respectively, represent the relative approximation error
$e_{U}$
for the propulsion speed of chemically and thermally active particles at different distances from the wall
$H/R$
. The inset shows the relative error for the drag
$e_{\hat{F}}$
as a function of
$H/R$
.
5 Conclusions
We theoretically examined the effect of confinement on the Marangoni propulsion of oblate spheroidal particles. We demonstrated that, depending on the geometry and degree of confinement, the surfing particles may propel forward in the higher surface tension direction or propel backward in the reverse direction. Knowing that the direction of propulsion is altered by a change in the surrounding boundary can be harnessed for designing smart surfing robots capable of sensing their environment.
In this work, we only considered particles with a
$90^{\circ }$
contact angle (i.e. half-submerged). However, the propulsion speed of a particle with a contact angle other than
$90^{\circ }$
can be estimated by calculating the speed of a half-submerged oblate spheroid with the same protrusion length to the liquid phase and the same contact-line radius. It has been shown that the fluid drags experienced by the two particles are closely comparable (see e.g. Stone & Masoud Reference Stone and Masoud2015).
Lastly, while our findings shed new light on the geometry–speed relation for Marangoni surfers, they also open up new questions such as how does the competition between the net surface tension force and the Marangoni flow generated by the surface tension gradients play itself out when the inertia is important, when the transport of the scalars is no longer dominated by diffusion or when the interface is curved? Answering these and other relevant questions will be the subject of future studies. It will also be interesting to see whether there exist examples of the reverse Marangoni propulsion in natural or biological systems.
Acknowledgements
We thank H. A. Stone for stimulating discussions. Partial support from the Nevada Advanced Autonomous Systems Innovation Center (NAASIC) is also acknowledged.