Capillary lubrication of a spherical particle near a fluid interface

Abstract

The lubricated motion of an object near a deformable boundary presents striking subtleties arising from the coupling between the elasticity of the boundary and lubricated flow, including but not limited to the emergence of a lift force acting on the object despite the zero Reynolds number. In this study, we characterize the hydrodynamic forces and torques felt by a sphere translating in close proximity to a fluid interface, separating the viscous medium of the sphere's motion from an infinitely more viscous medium. We employ lubrication theory and perform a perturbation analysis in capillary compliance. The dominant response of the interface owing to surface tension results in a long-ranged interface deformation, which leads to a modification of the forces and torques with respect to the rigid reference case, that we characterize in detail with scaling arguments and numerical integrations.

1. Introduction

The dynamics of objects moving in viscous fluids has been studied both theoretically and experimentally for a long time (Jeffery 1915; Lamb 1924; Collins 1955; Dean & O'Neill 1963; O'Neill 1964; Batchelor 1967). Confining the viscous flow between an object and a rigid surface modifies the forces felt by the object (Goldman, Cox & Brenner 1967; O'Neill & Stewartson 1967; Cooley & O'Neill 1969; Jeffrey & Onishi 1981). Such a modification is involved in vastly different phenomena ranging from the mechanics of joints (Hou et al. 1992; Hlaváček 1993), to the movement of cells in capillaries (Abkarian, Lartigue & Viallat 2002) and the dynamics of suspensions (Batchelor 1970, 1971, 1976, 1977; Happel & Brenner 1983).

Recent research has provided evidence of boundary elasticity further modifying the lubricated dynamics of an object (Leroy & Charlaix 2011; Leroy et al. 2012). Further standardization of the measurement process has led to the design of contactless probes for rheology (Garcia et al. 2016; Basoli et al. 2018). The coupling of boundary elasticity and lubrication flow, collectively termed soft lubrication, predicts the emergence of lift forces exerted on particles translating parallel to soft boundaries (Sekimoto & Leibler 1993; Beaucourt, Biben & Misbah 2004; Skotheim & Mahadevan 2005; Weekley, Waters & Jensen 2006; Urzay, Llewellyn Smith & Glover 2007; Snoeijer, Eggers & Venner 2013; Bouchet et al. 2015; Salez & Mahadevan 2015; Essink et al. 2021; Bertin et al. 2022; Bureau, Coupier & Salez 2023; Rallabandi 2024). Such lift forces are associated with the symmetry breaking arising out of the deformability. Since the latter is crucial to the generated force, the nature of the bounding wall has been further explored by examining the influence of slip (Rinehart et al. 2020), and viscoelasticity (Pandey et al. 2016; Kargar-Estahbanati & Rallabandi 2021). A reversal of the nature of the lift force from repulsive to attractive has also been predicted for viscoelastic settings (Hu, Meng & Doi 2023). Other studies have explored the complex modifications induced by including inertial effects (Clarke & Potnis 2011) and compressibility (Balmforth, Cawthorn & Craster 2010). On the experimental front, dedicated research has verified the presence of these lift forces on various substrates (Saintyves et al. 2016; Davies et al. 2018; Rallabandi et al. 2018; Vialar et al. 2019; Zhang et al. 2020).

In biology, where cells and tissues are extremely soft, and/or at small scales in soft matter, the interfacial capillary stress at the boundary dominates over bulk elasticity. By employing a classical Stokeslet-like response of the flow near a fluid interface, it has been shown that a rectified flow may be generated owing to the tension of the boundary (Aderogba & Blake 1978; Nezamipour & Najafi 2021). On the other hand, finite-size effects were addressed (Lee, Chadwick & Leal 1979; Lee & Leal 1980, 1982; Geller, Lee & Leal 1986) in the regime of a large gap between the object and the fluid interface, predicting counter-intuitive behaviours unique to capillarity. The results of these studies have been useful in analysing the movement of microorganisms near a fluid interface (Trouilloud et al. 2008; Lopez & Lauga 2014), as well as the formation of floating biofilms (Desai & Ardekani 2020). Recently, capillary-lubrication studies (Jha, Amarouchene & Salez 2023, 2024) have characterized the dynamics of an infinite cylinder near a fluid interface, highlighting the influence of the viscosity contrast and thickness ratio between the two fluid layers, and leading to large variability in the forces generated as opposed to elastic interfaces.

While previous research has highlighted the importance of understanding lubricated motion near a fluid interface, the characterization of the dynamics of a particle moving with multiple degrees of freedom near a fluid interface remains to be explored. In this article, we explore in detail the translational motion of a sphere moving in close proximity to an infinitely viscous but deformable sublayer. In the small-deformation limit, we calculate the forces and torques generated on the sphere during the motion. Due to a symmetry between translational and rotational motions in soft lubrication (Bertin et al. 2022), our work immediately generalizes to the case where rotation would be added. The remainder of the article is organized as follows. We start by describing the capillary-lubrication framework, before presenting the theoretical methodology for obtaining the different fields using perturbation analysis at small deformations of the fluid interface. We then discuss the implications of the interfacial deformation, and the competition between gravity and capillarity, on the forces and torques generated on the particle. Limiting expressions are derived for the capillary-dominated and gravity-dominated responses. While the former case is novel, the latter is reminiscent of a Winkler solid.

2. Capillary-lubrication theory

We consider a sphere of radius $a$ translating with a prescribed time-dependent horizontal velocity $\boldsymbol {u}=u(t)\boldsymbol {e}_x$ near a fluid interface, as shown in figure 1. The interface is characterized by its surface tension $\sigma$, and separates two incompressible Newtonian viscous liquids with dynamic shear viscosities $\eta$ and $\eta _{sl}$, as well as densities $\rho$ and $\rho _{sl}$ (with $\rho <\rho _{sl}$). The acceleration due to gravity is denoted by $g$. The gap profile between the sphere and the undeformed fluid interface is denoted by $h(r,t)$, which depends on the horizontal radial coordinate $r$ and time $t$. The $x$-direction oriented along $\boldsymbol {e}_x$ corresponds to the horizontal angular coordinate $\theta =0$. The time-dependent distance between the bottom of the sphere's surface and the undeformed interface is denoted by $d(t)$, with $d>0$. Its temporal derivative $\dot {d}(t)$ denotes the vertical velocity of the sphere. We focus on the case where the bottom layer is extremely viscous compared with the top layer, i.e. $\eta _{sl}\gg \eta$, and is infinitely thick compared with the gap between the sphere and the undeformed interface. These conditions allow us to assume that the tangential velocity at the fluid interface vanishes because the shear stress from the top layer is insufficient to generate a flow in the lower layer.

Figure 1. Schematic of the system. A sphere of radius $a$ immersed in a viscous fluid of viscosity $\eta$ and density $\rho$ moves near a fluid interface. The undeformed gap profile is noted $h(r,t)$, with $r$ the horizontal radial coordinate and $t$ the time. The origin of coordinates is located at the undeformed fluid interface ($z = 0$) in line with the centre of mass of the sphere ($r = 0$). The interface separates the top fluid from a secondary fluid of viscosity $\eta _{sl}$, with $\eta _{sl}\gg \eta$, and density $\rho _{sl}$ at the bottom, i.e. $\rho _{sl}>\rho$. The sphere has prescribed horizontal velocity $u$ and vertical velocity $\dot {d}$, where $d=h(0,t)$ denotes the instantaneous distance between the sphere bottom and the undeformed fluid interface. The interface deflection field is denoted as $\delta (r,t)$, the acceleration due to gravity is denoted as $g$ and the surface tension is denoted as $\sigma$.

We neglect fluid inertia and assume that the typical gap between the sphere and the interface is much smaller than the typical horizontal length scale, allowing us to invoke lubrication theory (Reynolds 1886; Oron, Davis & Bankoff 1997). Introducing the excess pressure field $p$, in the top layer, with respect to the hydrostatic contribution, and the horizontal velocity field $\boldsymbol {v}$ of the fluid in the gap, the incompressible Stokes equations thus read at leading lubrication order

(2.1)$$\begin{gather} \frac{\partial p}{\partial z} = 0, \end{gather}$$

(2.2)$$\begin{gather}\boldsymbol{\nabla} p= \eta\frac{\partial^2 \boldsymbol{v}}{\partial z^2}, \end{gather}$$

where $\boldsymbol {\nabla }$ denotes the gradient in the horizontal plane and where $z$ is the vertical coordinate. In the limit of a small gap, the thickness profile of the fluid between the sphere and the undeformed fluid interface can be approximated by its parabolic expansion, leading to

(2.3)\begin{equation} h(r,t) \simeq d(t)+\frac{r^{2}}{2 a}.\end{equation}

The no-slip boundary condition is imposed at both the sphere's surface and the fluid interface. Hence, in the frame of the moving sphere, at $z = h$, one has

(2.4)\begin{equation} \boldsymbol{v} = \boldsymbol{0}, \end{equation}

and at the interface, i.e. $z = \delta$, one has

(2.5)\begin{equation} \boldsymbol{v} ={-}\boldsymbol{u} ={-}u\boldsymbol{e}_x. \end{equation}

Given the boundary conditions above, we can easily write the horizontal velocity profile in the gap between the sphere and the interface as

(2.6)\begin{equation} \boldsymbol{v} = \frac{\boldsymbol{\nabla} p}{2\eta}(z-h)(z-\delta)+\boldsymbol{u}\frac{z-h}{h-\delta}.\end{equation}

The conservation of the fluid volume in the gap allows for the derivation of the Reynolds thin-film equation for the system, which reads

(2.7)\begin{equation} \frac{\partial}{\partial t}(h-\delta) = \boldsymbol{\nabla}\boldsymbol{\cdot} \left[\frac{\boldsymbol{\nabla} p}{12\eta}(h-\delta)^3+\frac{\boldsymbol{u}}{2}(h-\delta)\right].\end{equation}

The deflection of the fluid interface is controlled by the Laplace pressure jump. Thus, at $z = \delta$, one has

(2.8)\begin{equation} p \simeq\sigma \nabla^2\delta+g\delta(\rho-\rho_{sl}), \end{equation}

where we have assumed a small deformation of the interface and linearized the curvature.

Let us now non-dimensionalize the equations through

(2.9) \begin{equation} \left.\begin{gathered} d(t)= d^*D(T),\quad h(r,t) = d^*H(R,T),\quad \boldsymbol{r}=l\boldsymbol{R},\quad z =d^*Z,\\ t = \frac{l}{c}T,\quad \boldsymbol{v}(\boldsymbol{r},z,t)=c\boldsymbol{V}(\boldsymbol{R},Z,T),\quad \boldsymbol{u}(t) =cU(T)\boldsymbol{e}_x,\\ p(\boldsymbol{r},t) =\frac{\eta c l}{{d^*}^2}P(\boldsymbol{R},T),\quad \delta(\boldsymbol{r},t) = d^*\varDelta(\boldsymbol{R},T), \end{gathered}\right\} \end{equation}

where $c$ and $d^*$ represent characteristic in-plane velocity and gap-thickness scales, respectively, with $l = \sqrt {2ad^*}$ denoting the characteristic hydrodynamic radius.

At this point, we discuss the assumption of an infinitely viscous bottom layer and the resulting vanishing of the velocity, which has been previously addressed in the literature (Yiantsios & Davis 1990). Since the bottom fluid layer is considered to be much thicker than the upper one, the flow generated in the former evolves over a much larger length scale, typically given by $a$. As a consequence, by adapting the thin-bilayer-film normal stress balance (see (A4) in Bertin et al. 2021) to the current case of a thick bottom layer, we find that the dimensionless normal viscous stress at the fluid interface has a $(\eta _{sl}/\eta )(d^*/l)^4$ prefactor. We assume this numerical prefactor to be very small in our study, due to the vanishing of $d^*/l$, despite the large viscosity ratio – thus the non-inclusion of bottom-layer stress terms in (2.8).

In dimensionless terms, the undeformed thickness profile and Reynolds equation become

(2.10)\begin{equation} H(R,T) \simeq D(T)+R^{2}, \end{equation}

and

(2.11)\begin{equation} \frac{\partial}{\partial T}(H-\varDelta) = \boldsymbol{\nabla}\boldsymbol{\cdot} \left[\frac{\boldsymbol{\nabla} P}{12}(H-\varDelta)^3+\frac{\boldsymbol{U}}{2}(H-\varDelta)\right].\end{equation}

Besides, the deflection field $\varDelta$ is related to the excess pressure field $P$ by the dimensionless version of the Laplace equation, which reads

(2.12)\begin{equation} \nabla^2\varDelta-{Bo}\varDelta=\kappa P, \end{equation}

where ${Bo} = (l/l_{c})^2$ denotes the Bond number, $l_{c} = \sqrt {\sigma /[g(\rho _{sl}-\rho )]}$ denotes the capillary length and where we have introduced the capillary compliance of the interface, denoted by

(2.13)\begin{equation} \kappa = \frac{\eta cl^3}{\sigma {d^*}^3}. \end{equation}

3. Perturbation analysis

As in previous studies regarding soft lubrication (Sekimoto & Leibler 1993; Skotheim & Mahadevan 2005; Urzay et al. 2007; Snoeijer et al. 2013; Salez & Mahadevan 2015; Essink et al. 2021; Bertin et al. 2022), we perform an asymptotic expansion of the unknown fields at first order in dimensionless capillary compliance, $\kappa$, as

(3.1)$$\begin{gather} \varDelta \simeq \kappa \varDelta_{1}+O(\kappa^2), \end{gather}$$

(3.2)$$\begin{gather}P \simeq P_{0}+\kappa P_{1}+O(\kappa^2), \end{gather}$$

where $\kappa ^i\varDelta _i$ is the $i$th-order contribution to the deflection of the interface, and $\kappa ^iP_i$ is the $i$th-order contribution to the excess pressure field.

3.1. Zeroth-order pressure

At zeroth order in $\kappa$, (2.11) reduces to

(3.3)\begin{equation} \dot{D} = \boldsymbol{\nabla}\boldsymbol{\cdot} \left(\frac{\boldsymbol{\nabla} P_{0}}{12}H^3+\boldsymbol{U}\frac{H}{2}\right). \end{equation}

This equation is identical to the one for a perfectly rigid, flat and no-slip boundary. Since the zeroth-order pressure field results from linear terms in velocity, we decompose it azimuthally as

(3.4)\begin{equation} P_0(\boldsymbol{R},T) = P_{00}(R,T)+P_{01}(R,T)\cos\theta, \end{equation}

with $R$ and $\theta$ the horizontal polar coordinates of $\boldsymbol {R}$. Assuming a vanishing pressure field at large $R$, and that it must remain finite and single valued at $R=0$, the equations can be solved to give the zeroth-order components of the excess pressure field as (O'Neill & Stewartson 1967)

(3.5)$$\begin{gather} P_{00} ={-}\frac{3\dot{D}}{2(D+R^2)^2}, \end{gather}$$

(3.6)$$\begin{gather}P_{01} = \frac{6UR}{5(D+R^2)^2}. \end{gather}$$

3.2. Interface deflection

We now decompose the first-order deflection of the interface azimuthally as

(3.7)\begin{equation} \varDelta_1(\boldsymbol{R},T) = \varDelta_{10}(R,T)+\varDelta_{11}(R,T)\cos\theta. \end{equation}

Writing (2.12) at first order in $\kappa$, one has

(3.8)\begin{equation} \frac{1}{R}\frac{\partial}{\partial R}\left(R\frac{\partial \varDelta_{10}}{\partial R}\right)- {Bo}\varDelta_{10}=P_{00}, \end{equation}

and

(3.9)\begin{equation} \frac{1}{R}\frac{\partial}{\partial R}\left(R\frac{\partial \varDelta_{11}}{\partial R}\right)- \frac{\varDelta_{11}}{R^2}-{Bo}\varDelta_{11}=P_{01}, \end{equation}

with the boundary conditions $\varDelta _{1i}\rightarrow 0$ as $R\rightarrow \infty$, and finite and single-valued $\varDelta _{1i}$ at $R=0$. The solutions of the above equations can be written in the most general form as

(3.10)\begin{align} \varDelta_{10}(R) &={-}I_0(\sqrt{Bo}R)\int_R^\infty K_0(\sqrt{Bo}\xi)\xi P_{00}(\xi)\,\textrm{d}\kern0.7pt \xi \nonumber\\ &\quad -K_0(\sqrt{Bo}R)\int_0^R I_0(\sqrt{Bo}\xi)\xi P_{00}(\xi)\,\textrm{d}\kern0.7pt \xi, \end{align}

(3.11)\begin{align} \varDelta_{11}(R) &={-}I_1(\sqrt{Bo}R)\int_R^\infty K_1(\sqrt{Bo}\xi)\xi P_{01}(\xi)\,\textrm{d}\kern0.7pt \xi \nonumber\\ &\quad -K_1(\sqrt{Bo}R)\int_0^R I_1(\sqrt{Bo}\xi)\xi P_{01}(\xi)\,\textrm{d}\kern0.7pt \xi, \end{align}

where $I_j$ and $K_j$ denote the $j$th-order modified Bessel functions of the first and second kinds, respectively.

To understand the parametric influence of the Bond number ${Bo}$, we explore the behaviours of the interface deflection for both small and large ${Bo}$ values. For vanishing ${Bo}$, the anisotropic deflection component $\varDelta _{11}$ reaches a limiting behaviour given by the following expression:

(3.12)\begin{equation} \varDelta_{11} \simeq{-}\frac{3U}{10}\frac{\ln\left(1+\dfrac{R^2}{D}\right)}{R}. \end{equation}

In contrast, the isotropic deflection component $\varDelta _{10}$ does not show any limiting behaviour at vanishing ${Bo}$, and a reduction in ${Bo}$ leads to an unbounded increase in $\varDelta _{10}$. To understand the behaviour of $\varDelta _{10}$ as ${Bo}\rightarrow 0$, we take an asymptotic approach previously used in problems relating to capillary deformations (James 1974; Lo 1983; Dupré de Baubigny et al. 2015). The vanishing of ${Bo}$ allows for a scale separation in the radial coordinate $R$ into: (i) an inner region controlled by surface tension; and (ii) an outer region, where gravity is present. In the inner region ($R\ll 1/\sqrt {Bo}$), the interface deformation denoted by $\varDelta _{10}^{i}$ satisfies

(3.13)\begin{equation} \frac{1}{R}\frac{\partial}{\partial R}\left(R\frac{\partial \varDelta_{10}^{i}}{\partial R}\right)=P_{00}. \end{equation}

The general solution of the above equation reads

(3.14)\begin{equation} \varDelta_{10}^{i} = \mathcal{A}-\frac{3\dot{D}}{8D}\ln(D+R^2), \end{equation}

where $\mathcal {A}$ is an integration constant. The far-field behaviour of the inner solution reads

(3.15)\begin{equation} \varDelta_{10}^{i} \sim \mathcal{A}-\frac{3\dot{D}}{4D}\ln(R). \end{equation}

In the outer region ($R \gg 1/\sqrt {Bo}$), gravity matters but there is no lubrication pressure, which leads to the following governing equation for the interface deformation denoted by $\varDelta _{10}^{o}$:

(3.16)\begin{equation} \frac{1}{R}\frac{\partial}{\partial R}\left(R\frac{\partial \varDelta_{10}^{o}}{\partial R}\right)- {Bo}\varDelta_{10}^{o}=0. \end{equation}

The latter equation is solved, along with the condition that the deflection vanishes at infinity, leading to

(3.17)\begin{equation} \varDelta_{10}^{o} = \mathcal{B}K_0(\sqrt{Bo}R), \end{equation}

where $\mathcal {B}$ is an integration constant. The near-field behaviour of the outer solution reads

(3.18)\begin{equation} \varDelta_{10}^{o} \sim{-}\mathcal{B}\left(\gamma-\ln 2+\frac{\ln {Bo}}{2}+\ln R\right), \end{equation}

where $\gamma$ is Euler's constant. Matching (3.15) with (3.18) leads to

(3.19)$$\begin{gather} \mathcal{A} ={-}\frac{3\dot{D}}{4D}\left(\gamma-\ln 2+\frac{\ln {Bo}}{2}\right), \end{gather}$$

(3.20)$$\begin{gather}\mathcal{B} = \frac{3\dot{D}}{4D}. \end{gather}$$

Substituting these constants in (3.14) and (3.17), we find the matched asymptotes of the interface deflection at small ${Bo}$. The interface deflection can then be approximated by the matched cross-over expression

(3.21)\begin{equation} \varDelta_{10}|_{{Bo}\rightarrow0} \approx \frac{3\dot{D}}{4D} \left[K_0(\sqrt{Bo}R)+\frac12\ln\left(\frac{R^2}{D+R^2}\right)\right].\end{equation}

The interface deflection is shown in figure 2 at a fixed low value of ${Bo}$. The cross-over expression at small ${Bo}$ described above matches the exact one calculated using (3.10), with improving precision as ${Bo}$ reduces. Figure 2 also shows the inner and outer solutions, which diverge in the far and near fields, respectively.

Figure 2. Isotropic component $\varDelta _{10}$ of the amplitude of the first-order interface deflection as a function of the radial coordinate $R$ (solid black line), as calculated from (3.10), for ${Bo} = 0.01$, $D = 1$ and $\dot {D} = 1$. For comparison, the inner solution (red solid line), the outer solution (blue solid line) and the matched cross-over expression (symbols), from (3.14), (3.17) and (3.21), respectively, are also shown.

The other interesting limit arises when ${Bo}\rightarrow \infty$, leading to the curvature-related terms in (3.8)–(3.9) dropping out, giving us

(3.22)$$\begin{gather} \varDelta_{10} ={-}\frac{P_{00}}{Bo} = \frac{3\dot{D}}{2{Bo}(D+R^2)^2}, \end{gather}$$

(3.23)$$\begin{gather}\varDelta_{11} ={-}\frac{P_{01}}{Bo} ={-}\frac{6UR}{5{Bo}(D+R^2)^2}. \end{gather}$$

As a consequence, the deflection is directly proportional to the pressure applied, with $\kappa /{Bo} = \eta c l/[g{d^*}^3(\rho _{sl}-\rho )]$ as an effective compliance. From the latter, we see that in the limit of large ${Bo}$ the surface tension no longer controls the deformation. This large-${Bo}$ response is akin to the Winkler response for thin compressible elastic materials, which has been studied previously (Skotheim & Mahadevan 2005; Urzay et al. 2007; Salez & Mahadevan 2015; Chandler & Vella 2020; Bertin et al. 2022).

Apart from these two extremes, further exploration of the influence of ${Bo}$ can be done using (3.10) and (3.11). Figure 3 shows the isotropic and anisotropic profiles of the amplitude of the first-order interface deflection, for different values of ${Bo}$. We see that capillarity leads to a long-ranged interface deflection whose range and magnitude both decrease with increasing ${Bo}$.

Figure 3. Isotropic component $\varDelta _{10}$ (a) and anisotropic component $\varDelta _{11}$ (b) of the amplitude of the first-order interface deflection as a function of the radial coordinate $R$, as calculated from (3.10) and (3.11), for $D = 1$, various ${Bo}$ as indicated, and for either $\dot {D} = 1$ (a) or $U=1$ (b). The black solid line denotes the limiting profile for ${Bo} = 0$ in the anisotropic case, as given in (3.12).

3.3. First-order pressure

At first order in $\kappa$, the pressure field involves in particular the squared velocity, and can thus be decomposed as

(3.24)\begin{align} P_1(\boldsymbol{R},T,{Bo},D) = P_{10}(R,T,{Bo},D)+P_{11}(R,T,{Bo},D)\cos \theta+P_{12}(R,T,{Bo},D)\cos2\theta. \end{align}

The governing equations for the components $P_{1i}$ of the first-order magnitude $P_1$ of the excess pressure field can be derived by considering (2.11) at first order in $\kappa$. Since $P_{12}$ is not needed to compute the forces and torques, we restrict ourselves to the two following equations:

(3.25)\begin{align} \frac{1}{R}\frac{\partial}{\partial R}\left(RH^3\frac{\partial P_{10}}{\partial R}\right) &= \frac{1}{R}\frac{\partial}{\partial R}\left[3RH^2\left(\varDelta_{10} \frac{\partial P_{00}}{\partial R}+\frac{\varDelta_{11}}{2} \frac{\partial P_{01}}{\partial R}\right)\right] \nonumber\\ &\quad +3U\frac{\partial\varDelta_{11}}{\partial R} +3U \frac{\varDelta_{11}}{R}-12\frac{\partial \varDelta_{10}}{\partial T}, \end{align}

(3.26)\begin{align} \frac{1}{R}\frac{\partial}{\partial R}\left(RH^3\frac{\partial P_{11}}{\partial R}\right)- \frac{H^3}{R^2}P_{11} &= \frac{1}{R}\frac{\partial}{\partial R}\left[3RH^2\left(\varDelta_{10} \frac{\partial P_{01}}{\partial R}+\varDelta_{11}\frac{\partial P_{00}}{\partial R}\right)\right] \nonumber\\ &\quad -\frac{3H^2\varDelta_{10}P_{01}}{R^2}+6U\frac{\partial\varDelta_{10}}{\partial R}-12\frac{\partial \varDelta_{11}}{\partial T}. \end{align}

Using the linearity of the equations above, the components of the pressure field can be expressed as follows:

(3.27)$$\begin{gather} P_{10} = \frac{U^2}{D^2}\phi_{U^2}+\frac{\dot{D}^2}{D^3}\phi_{\dot{D}^2}+ \frac{\ddot{D}}{D^2}\phi_{\ddot{D}}, \end{gather}$$

(3.28)$$\begin{gather}P_{11} = \frac{\dot{U}}{D}\phi_{\dot{U}}+\frac{U\dot{D}}{D^2}\phi_{U\dot{D}}, \end{gather}$$

where the $\phi _k$ represent the auxiliary functions for the corresponding second-order forcing parameters $k$, such as $U^2$ etc., which all vary with distance $R$ and depend upon the values of $D$ and ${Bo}$. These functions can then be evaluated by numerically solving (3.25)–(3.26), together with the boundary conditions that the pressure vanishes in the far field and remains finite and single valued at the origin. The results are shown in figures 4 and 5. For all components, as ${Bo}$ increases, the magnitudes of the auxiliary functions decrease due to the decreasing interface deflection. Interestingly, even though the deflection studied above is long ranged, the pressure field decays quite rapidly.

Figure 4. Auxiliary functions $\phi _k$ (see (3.27)) of the isotropic component $P_{10}$ of the first-order magnitude $P_1$ of the excess pressure field, as functions of the radial coordinate $R$, for various Bond numbers ${Bo}$, as obtained by numerically solving (3.25) with $D = 1$.

Figure 5. Auxiliary functions $\phi _k$ (see (3.28)) of the anisotropic component $P_{11}$ of the first-order magnitude $P_1$ of the excess pressure field, as functions of the radial coordinate $R$, for various Bond numbers ${Bo}$, as obtained by numerically solving (3.26) with $D = 1$.

4. Capillary-lubrication forces and torques

Since the zeroth-order forces and torques acting on the sphere are identical to the known ones for the motion near a rigid, flat and no-slip boundary, we focus here on the first-order forces and torques acting on the sphere, and resulting from the interface deflection. These can be calculated from the stress tensor $\boldsymbol {\varSigma }$, that reads in the lubrication approximation: $\boldsymbol {\varSigma } \simeq -p\boldsymbol {I}+\eta \boldsymbol {e}_z\partial _z\boldsymbol {v}$, where $\boldsymbol {I}$ denotes the identity tensor. In dimensional units, the first-order vertical force acting on the sphere can be evaluated as

(4.1)\begin{equation} f_{z,1} \simeq \frac{2^{5/2}{\rm \pi}\kappa\eta ca^{3/2}}{{d^*}^{1/2}}\int_0^{\infty} P_{10}(R) R\,\textrm{d}R, \end{equation}

which can be decomposed using the auxiliary functions calculated before as

(4.2)\begin{equation} f_{z,1} \simeq \alpha_{U^2}\frac{\eta^2 u^2 a^3}{\sigma d^2}-\alpha_{\dot{D}^2} \frac{\eta^2 \dot{d}^2 a^4}{\sigma d^3}+\alpha_{\ddot{D}}\frac{\eta^2 \ddot{d} a^4}{\sigma d^3}, \end{equation}

where the $\alpha _k$ (with $k$ indicating here the forcing source, such as $U^2$) are the prefactors of the respective scalings. These prefactors are plotted in figure 6 as functions of ${Bo}$. An important difference arises between the various prefactors at small ${Bo}$. Indeed the prefactor $\alpha _{U^2}$ reaches a plateau, whereas $\alpha _{\ddot {D}}$ and $\alpha _{\dot {D}^2}$ increase logarithmically with decreasing ${Bo}$. These asymptotic behaviours at small ${Bo}$ can be calculated using the Lorentz reciprocal theorem (Daddi-Moussa-Ider et al. 2018; Rallabandi et al. 2018; Masoud & Stone 2019), by invoking as well (3.12) and (3.21), as detailed in Appendix A and summarized in table 1. They are in agreement with the numerical solutions, as shown in figure 6. Besides, as ${Bo}$ increases, all the prefactors decrease, and this decrease becomes inversely proportional to ${Bo}$ in the large-${Bo}$ limit, highlighting the transition to the Winkler-like regime. The corresponding asymptotic expressions have been calculated previously (Bertin et al. 2022), are summarized in table 1, and are in agreement with the numerical solutions, as shown in figure 6.

Figure 6. Prefactors $\alpha _k$ of the scalings of the vertical-force terms, defined in (4.2), as functions of Bond number ${Bo}$. The solid and dashed lines correspond to the small- and large-${Bo}$ behaviours (see table 1), respectively.

Table 1. Asymptotic behaviours of the scaling prefactors $\alpha _k$ and $\beta _k$ at small (see Appendix A) and large (see Bertin et al. 2022) Bond numbers ${Bo}$, where $\gamma$ denotes Euler's constant.

Similarly, the horizontal force acting on the sphere is given by the expression (Bertin et al. 2022)

(4.3)\begin{align} f_{x,1} &\simeq 2{\rm \pi}\eta c a\kappa\int_0^\infty \left[{-}2RP_{11}-\frac{H}{2}\left(\partial_R P_{11}+\frac{P_{11}}{R}\right)+\frac{\varDelta_{11}}{2}\partial_R P_{00} \right. \nonumber\\ &\quad \left. +\frac{\varDelta_{10}}{2}\left(\partial_RP_{01}+\frac{P_{01}}{R}\right)-2 \frac{U\varDelta_{10}}{H^2}\right]R\,\textrm{d}R, \end{align}

which can be decomposed into

(4.4)\begin{equation} f_{x,1} \simeq{-}\beta_{U\dot{D}}\frac{\eta^2 u \dot{d}a^3}{\sigma d^2}+\beta_{\dot{U}} \frac{\eta^2 \dot{u}a^3}{\sigma d},\end{equation}

where the $\beta _k$ (with $k$ indicating here the forcing source, such as $\dot {U}$) are the prefactors of the respective scalings. These prefactors are plotted in figure 7 as functions of ${Bo}$. An important difference arises once again between the two prefactors at small ${Bo}$. Indeed, the prefactor $\beta _{\dot {U}}$ reaches a plateau, whereas $\beta _{U\dot {D}}$ increases logarithmically with decreasing ${Bo}$. These asymptotic behaviours at small ${Bo}$ can be once again calculated using the Lorentz reciprocal theorem, as detailed in Appendix A and summarized in table 1. They are in agreement with the numerical solutions, as shown in figure 7. Besides, as ${Bo}$ increases, all the prefactors decrease, and this decrease becomes inversely proportional to ${Bo}$ in the large-${Bo}$ limit, highlighting once again the transition to the Winkler-like regime. The corresponding asymptotic expressions have been calculated previously (Bertin et al. 2022), are summarized in table 1 and are in agreement with the numerical solutions, as shown in figure 7.

Figure 7. Prefactors $\beta _k$ of the scalings of the horizontal-force terms, defined in (4.4), as functions of Bond number ${Bo}$. The solid and dashed lines correspond to the small- and large-${Bo}$ behaviours (see table 1), respectively.

As shown in previous studies (Salez & Mahadevan 2015; Bertin et al. 2022), the contributions to the first-order torque felt by the sphere along the $y$-axis have the same numerical prefactors as those for the first-order horizontal force, with the inclusion of a supplementary length-scale factor $a$. Hence, the first-order torque exerted on the sphere is given by

(4.5)\begin{equation} \tau_{y,1} \simeq \beta_{U\dot{D}}\frac{\eta^2 u \dot{d}a^4}{\sigma d^2}-\beta_{\dot{U}} \frac{\eta^2 \dot{u}a^4}{\sigma d}.\end{equation}

We conclude this whole section with an important remark. The cross-over from the capillary-dominated to the Winkler-like regime occurs at ${Bo} \sim 1$. This occurs when the hydrodynamic radius is comparable to the capillary length. Since typical capillary lengths are of the order of $\sim$1 mm, and given the lubrication conditions, such a cross-over can only be felt with spheres of millimetric radii and above.

5. Conclusion

Using capillary-lubrication theory, scaling arguments and numerical integrations, we explored the asymptotic forces and torques generated on a sphere in translational motion within a viscous fluid, in close proximity to a deformable capillary interface with another, infinitely viscous fluid on the other side. Due to a symmetry between translational and rotational motions in soft lubrication (Bertin et al. 2022), our work immediately generalizes to the case where rotation would be added. Specifically, by employing a perturbation analysis in the limit of small deformation of the fluid interface, we calculated the pressure fields decomposed into their various contributions from different degrees of freedom of the sphere. We investigated in particular the effects of gravity, which not only change the scaling laws of the forces and torques, but also show a Winkler-like elastic response at large Bond numbers. Altogether, our results allow us to quantify and possibly control soft-lubricated motion near tensile interfaces, which are ubiquitous in soft matter and biological physics. For example, the visco-capillary lift force revealed among others by our analysis, could contribute to the rationalization of the swimming of active organisms at the air–water interface (Trouilloud et al. 2008), the self-propulsion of levitating Leidenfrost droplets (Gauthier et al. 2019) for which a detailed mechanism is lacking to date or the thermal motion of oil droplets near rigid walls for which intriguing transient non-conservative forces have been reported (Fares et al. 2024).