We investigate the flow of viscous interfaces carrying an insoluble surface active material, using numerical methods to shed light on the complex interplay between Marangoni stresses, compressibility, and surface shear and dilatational viscosities. We find quantitative relations between the drag on a particle and interfacial properties as they are required in microrheology, i.e., going beyond the asymptotic limits. To this end, we move a spherical particle probe at constant tangential velocity, symmetrically immersed at either the incompressible or compressible interface, in the presence and absence of surfactants, for a wide range of system parameters. A full three-dimensional finite element calculation is used to reveal the intimate coupling between the bulk and interfacial flows and the subtle effects of the different physical effects on the mixed-type velocity field that affects the drag coefficient, both in the bulk and at the interface. For an inviscid interface, the directed motion of the particle leads to a gradient in the concentration of the surface active species, which in turn drives a Marangoni flow in the opposite direction, giving rise to a force exerted on the particle. We show that the drag coefficient at incompressible interfaces is independent of the origin of the incompressibility (dilatational viscosity, Marangoni effects or a combination of both) and that its higher value can not only be related to the Marangoni effects, as suggested earlier. In confined flows, we show how the interface shear viscosity suppresses the vortex at the interface, generates a uniform flow, and consequently increases the interface compressibility and the Marangoni force on the particle. We mention available experimental data and provide analytical approximations for the drag coefficient that can be used to extract surface viscosities.

## I. INTRODUCTION

The rheology of complex fluid-fluid interfaces is a relevant and rich physical problem.^{1,2} For example, the seminal theoretical work of Saffman and Delbrück^{3} on the Brownian motion in biological membranes revealed how the momentum diffusion at the interface and in the bulk are strongly coupled and, for non-compressible interfaces, the ratio between bulk and surface viscosities determines the range of the hydrodynamic effects at the interface. However, when Brownian motion is used to more generally deduce interfacial rheological properties beyond lipid membranes, values obtained from microrheology have been found to be very different from those with macroscopic interfacial shear rheometers.^{4} Several factors have been identified to contribute to the drag of particles close to and at an interface, and possible sources of discrepancies have been identified, with some being related to dissipation modes of the surface or to electrokinetic effects^{5} and others focused on the non-correct analysis of the hydrodynamic conditions of the measurements,^{6} although for direct shear rheological techniques such hydrodynamic analysis has led to clear operating windows.^{7,8} One aspect which merits further investigation is the role of interfacial compressibility and the possible influence of Marangoni flows on the interfacial and bulk velocity fields. First, there is clear experimental evidence from studying needles of different aspect ratios that compression might be of relevance in interpreting microrheological experiments.^{9} Second, the analysis of different types of interfacial rheometers, which use mixed interfacial flow fields for measuring the surface rheological properties in mixed shear-dilation/compression^{10} or shear-extension-dilation/compression,^{11} has already revealed a complex and subtle interplay between the Marangoni stresses, bulk flow effects, and the interfacial stresses stemming from the desired rheological material functions. In the present work, we investigate to what extent this analysis holds true for the flow field caused by a particle dragged through the interface (the disturbance velocity field), and resolving the effects of the various coupled contributions to the forces and concentration fields, as such a translating probe produces a mixed-type flow field.^{12} In particular, we are going to investigate the role of interfacial compressibility and the interplay between the different interfacial deformation modes (shear and local dilatation/compression).

The drag of an object at an interface has already received attention previously, and typically, two effects have been studied, i.e., the effect of coupling of the interfacial momentum diffusion to that in the bulk and the effect of Marangoni contributions. The original work of Saffman and Delbrück^{3} examined a cylindrical probe translating at a membrane in the limit of large membrane viscosity and an incompressible interfacial layer. In this limit, a dependence of drag on the particle size was found to be essentially absent (only a weak logarithmic effect had been predicted), and the range of the hydrodynamic forces was identified as the Saffman–Delbrück length. Later, Hughes *et al.*^{13} obtained the exact solution to the equations of motion for any combination of membrane and adjacent liquid viscosities. In both these works, the domain of the bulk flow was considered infinite. The influence of the bulk layer thickness was studied by Evans and Sackmann.^{14} They imposed a proportionality between the velocity and the shear stress exerted by the subphase on the membrane. This work has been generalized by Stone and Ajdari^{15} to any value of subphase depth using numerical calculations. Barentin *et al.*^{16,17} used the same geometry as Evans and Sackmann^{14} and studied both cases of an incompressible Langmuir (insoluble) monolayer and a Gibbs (soluble) monolayer. They exploited the lubrication approximation for shallow films and obtained analytical expressions for the shear stress and the drag force. In all studies mentioned so far, the particle was cylindrical and non-protruding.

For a protruding spherical particle trapped within a very thin, Newtonian fluid-fluid interface, the drag has been calculated only numerically. While Danov *et al.*^{18,19} considered a compressible interface characterized by a Boussinesq–Scriven model with shear and dilational viscosities, Fischer^{20} argued that due to the short time scales on which a sustained gradient in surface tension gives rise to Marangoni contributions to the forces, the monolayer can under certain conditions be considered incompressible, as the presence of a surfactant strongly suppresses any motion at the surface that compresses or expands the interface. Manikantan and Squires^{21} provided an exhaustive, excellent discussion of these phenomena. In the study of monolayers of fatty acids and phospholipids, the surface shear viscosity calculated following Danov's approach was shown to underestimate the experimental value.^{22}

Fischer *et al.*^{23} then solved the Stokes flow equations for a sphere immersed in a viscous interface (monolayer) upon incorporating the Marangoni effect into the treatment, by solving the equations for an incompressible interface. Fischer showed that the Marangoni effect contributes a significant part of the total drag at the limit of vanishing surface compressibility. However, not all systems will adhere to the conditions required for a vanishing surface compressibility, and how this will come into play is ill-understood, in particular when surface viscosities are significant as well. A recent approach by Elfring *et al.*^{12} to characterize the latter used the lubrication approximation for a very shallow subphase and considered a non-protruding disk-shaped probe at a viscous, compressible interface that generates Marangoni flows. Elfring *et al.* assumed the presence of a soluble surfactant that equilibrates on a finite timescale and did not consider transient effects. While the lubrication approximation allows for analytical solution through a perturbation expansion, using a shallow subphase is not ideal for probing rheological properties of the interface, because in this regime, bulk stresses dominate over surface stresses.^{12}

The present work can be considered as an extension of the work of Elfring *et al.* to a protruding particle and a deep subphase. We study the transient Marangoni flow and the interface compressibility of a Newtonian interface endowed with both shear and dilatational viscosity. The interface is covered by an insoluble surface active component (surfactant, protein, fatty acid, lipid, or polymeric substance) and a large spherical particle, which is set to translate at the interface. The interfacial contribution to the drag on the particle is a combination of the Marangoni effect and the interfacial stresses arising from the interfacial strain fields and relevant interfacial viscosities. This interfacial disturbance velocity field depends on both the balance between the Marangoni flows and the shear and local dilational/compressional components and how the momentum sink of the bulk influences this. To the best of our knowledge, this is the first study on a protruding three-dimensional particle at compressible interfaces with Marangoni flows. Dimova *et al.*^{24} also considered concentration variations with a spherical particle, but did not couple the interface to a bulk fluid and did not determine the effects of Marangoni flow.

We begin by recalling the hydrodynamic equations and formulate the boundary and initial conditions for the problem at hand in Sec. II. We then demonstrate how we calculate the drag on the particle and introduce dimensionless quantities to come up with reduced equations and a number of dimensionless numbers characterizing all parameters of the particle, the interface, and the bulk phases. The numerical implementation is described in Sec. III. In Sec. IV, we present results for various regimes, the surfactant dominated regime (Sec. IV A), a shear-free interface without (Sec. IV B), and with Marangoni flow (Sec. IV C). In Sec. IV D, we focus on the flow field of regime of finite shear viscosity, as well as a reduced parameter space that captures the situation encountered in real surfactant systems. The effects of interface viscosites on the drag coefficient are discussed in Sec. IV E. Conclusions are provided in Sec. V, with a caveat regarding the use of these drag flows for determining shear viscosities or rheological properties.

## II. HYDRODYNAMIC EQUATIONS

Consider the rather general situation of a spherical particle symmetrically immersed at a planar interface, with its equator residing in the interface plane, and translating with a constant velocity $U=Uex$ (Fig. 1). The two fluids in domains $\Omega 1$ and $\Omega 2$ are considered incompressible and Newtonian, i.e., modeled by

where ** r** denotes position, $\pi i$ is the stress tensor in domain

*i*,

*p*is the pressure field,

**is the unity tensor, $\tau i=\eta i[\u2207u+(\u2207u)T]$ is the viscous stress tensor,**

*I***is the fluid velocity, and**

*u**η*is the viscosity of fluid $i\u2208{1,2}$. All fields appearing in our equations, such as

_{i}**and $\tau i$, are functions of spatial coordinate**

*u***and time**

*r**t*, and we skip these arguments for better readability. A similar decomposition can also be written for the interface stress tensor

^{1,25}

where $\tau s$ is the extra surface stress tensor, $Is=I\u2212nn$ is the surface or tangential projection tensor^{25,26} with surface normal ** n,** and

*γ*is the surface tension of the interface. For an isothermal system, the surface tension is solely a function of the surfactant concentration Γ, i.e., $\gamma =\gamma (\Gamma )$. Its functional form remains to be specified (Sec. IV). The extra surface stress tensor is described by the Boussinesq–Scriven constitutive law

^{27,28}

where *η ^{s}* and

*κ*are the shear and dilatational viscosities of the interface, respectively. The surface viscosities can be made functions of the surface pressure $\Pi s(\Gamma )=\gamma 0\u2212\gamma (\Gamma )$, where

^{s}*γ*

_{0}is the surface tension of a clean interface (Γ = 0). It introduces additional complexities, such as irreversible particle motion.

^{29}However, for simplicity in the current analysis, they are here assumed to be system parameters, and an explicit dependence on Γ is not included.

^{28,30}In most practical cases, one needs some surface active species at the interface to get significant extra (viscous) stresses,

^{31,32}but surface viscosities for clean interfaces have also been reported.

^{33}The surface rate of deformation tensor $Ds$ is defined as

where $us$ is the velocity ** u** evaluated at the interface, and $\u2207s=Is\xb7\u2207$ is the surface gradient operator.

^{26}Note that this is a pseudo-definition of $\u2207s$, as $\u2207$ does not exist for a surface field. Given the small length scale of the particle involved, the dynamics of the incompressible bulk fluids can be modeled with the creeping-flow (or Stokes) equations, which reduce the momentum balance to $\u2207\xb7\pi i=0$. With the help of Eq. (1), this becomes

where the vanishing divergence of ** u** expresses the balance of mass assuming a constant density within the domains. The velocity and pressure fields thus receive their time-dependency through Γ. The evolution of the surfactant concentration Γ is governed by the unsteady surface convection-diffusion (SCD) equation

^{26,34–36}

where *D _{s}* is the surface diffusivity of the surfactant and where we assumed a planar interface. Equations (5) and (6) are the governing equations for

**and Γ as a function of time and**

*u***. They require proper boundary and initial conditions to be solved.**

*r*In this study, it is assumed that the interface remains planar, which is a valid assumption for the typical particle sizes and interfacial tensions used in experiments.^{37} Moreover, we solve the hydrodynamic problem only for the translational motion of the spherical colloids along the interface, i.e., under the assumption of contact line pinning, which prevents the particle from rotating. It has been shown experimentally that the presence of surface roughness leads to pinning of the contact line.^{38} The pinned contact line is also used for modeling single^{37} and pair^{39} particles at the interfaces. Due to pinning of the contact line, an extremely small deformation of the interface is needed to balance the torque acting on the particle, which justifies our approach of a non-rotating particle at a planar interface.

### A. Boundary and initial conditions

For the surfactant concentration Γ governed by Eq. (6), we use the following initial and boundary conditions:

where *S _{I}* is the fluid-fluid interface, $\u2202Sp$ is the circular particle-interface boundary, $\u2202Sb$ is the square box boundary, and $np$ is the unit vector normal to the particle surface and tangential to the interface. We thus assume an initially undisturbed, uniform Γ profile at the fluid-fluid interface, which remains at the initial value at all times at the box boundary. The third boundary condition means no surfactant flux from the particle interface boundary.

The boundary conditions for the velocity field ** u** are composed of boundary conditions on the surface

*S*of the particle, at the fluid-fluid interface

_{p}*S*, and far away from the particle on the simulation box surface

_{I}*S*. In this study, we assume that there is no-slip at the particle-fluid interface as well as at the interface between two Newtonian fluids, while more generally, there can be relative motion between surfaces in contact.

_{b}^{25,40,41}However, the slip boundary condition at the interface is usually only employed for macromolecular, non-Newtonian fluids, and its necessity for Newtonian fluids remains unclear.

^{40}We assume that the particle translates with velocity $U=r\u0307p$, where $rp$ is the location of the center-point of the particle, and the particle does not rotate.

^{42}Applying no-slip on the particle surface

*S*and vanishing velocity on the simulation box surface

_{p}*S*yields

_{b}Boundary conditions for ** u** at the fluid-fluid interface are defined by jump equations. In the absence of mass transfer (phase change) across the interface, there are the mass jump balance

^{25}

and the momentum jump balance

where $n=ey$ is interface surface normal pointing to region $\Omega 2$ (Fig. 1), $ui$ and *p _{i}* denote the velocity and pressure of fluid

*i*at $r\u2208SI$, and $\u2207s\xb7n$ is the interface mean curvature. The subscripts 1 and 2 indicate the value of quantity in the bulk phases extrapolated to the interface. The momentum jump balance can be split into normal,

and tangential components,

where ** t** is a unit tangent vector residing in the

*x*-

*z*-plane (Fig. 1).

In our setup, the interface remains at *y* = 0, coinciding with the *y*-coordinate of center of the sphere, the normal component of $us$ therefore vanishes, and the momentum jump balance in normal direction is replaced by $us\xb7n=0$. For the tangential component of $us$, we assume a no-slip boundary conditions at the interface, which implies continuity of the tangential component of the velocity across the interface,

Finally, we note that the tangential momentum jump balance can be written in in terms of the surface pressure^{25}

The Gibbs–Marangoni modulus $K\pi =\Gamma \u2202\Pi s/\u2202\Gamma $ allows one to relate the gradient in surface pressure to the gradient in concentration,^{21,26}

and to then rewrite the tangent component of the stress jump equation in the following final form:

### B. Drag on the particle

The general relation between the drag on the particle ** F** at the interface translating with a constant velocity

**is complex and depends on the particle geometry, bulk and interface rheological properties, the interfacial equation of state, and the transport properties of the surfactant. The drag coefficient on a spherical particle with radius**

*U**R*in an unbounded Newtonian fluid translating with a constant velocity is the Stokes' drag coefficient $f=|F|/|U|=6\pi \eta R$.

^{15}For a spherical particle symmetrically embedded in a clean, inviscid interface between a Newtonian and an inviscid fluid, the drag is exactly half the Stokes' drag, i.e., $f=3\pi \eta R$.

^{43}

The drag force on the spherical particle embedded at the interface is the sum of interface and bulk contributions,

where $\u2202Sp$ is the perimeter of the particle at the interface and $np$ is the unit normal vector to the surface of the particle (Fig. 1). The first integral on the right hand side of (17) is the bulk force $Fb$. The second integral is the interface force $Fi$ and can be decomposed into two contributions, i.e.,

The first integral on the right hand side of Eq. (18) is the elastic or Marangoni part of the drag force $FM$ due to the non-uniform distribution of the surfactant in the contact line.^{24} The second integral is the interface viscous force $Fs$ due to the extra interface stress tensor. Hence, the total drag force (17) on the particle is

Because our particle is non-rotating, a torque will act on the particle, which can be evaluated via suitable adapted integrals involving the relative position $x=r\u2212rp$, where $rp$ denotes the position of the particle's center. As explained earlier, it is assumed that this torque is balanced by an extremely small deformation of the interface.

### C. Dimensionless quantities and equations

For the problem at hand, we can use the constant velocity $U=|U|$, the radius *R* of the spherical particle, the viscosity *η*_{1} of fluid 1, and the initial surfactant concentration $\Gamma 0$ to introduce reduced units, and to come up with a number of dimensionless parameters. Dimensionless variables, marked by asterisk (only here), are therefore defined uniquely in terms of their dimensional counterparts, such as

From now on, and for the rest of this work, variables are meant to be dimensionless, and we omit the asterisk for brevity. Using them in the above sections, Eq. (5) is replaced by

where $\lambda =\eta 2/\eta 1$ is the ratio between the two viscosities. The normal and tangential parts of the jump balance, Eqs. (11) and (16) remain formally unaltered, i.e.,

with $\tau 1=2D$ and $\tau 2=2\lambda D$. To close the system of equations, and to fully specify the dimensionless quantities, an interfacial equation of state must be chosen, which will be done in Sec. IV. Because we are going to assume (i) a linearized equation of state there, $K\pi $ will be replaced by $Ma\u2009\Gamma $, where Ma is a dimensionless Marangoni number, (ii) a flat interface, the last term in the first line of Eq. (22) will disappear, and (iii) an air–liquid interface, *λ* = 0.

The dimensionless extra surface stress tensor (3)

is written in terms of the Bousinessq number $Bq1=\eta s/R\eta 1$, the ratio of surface shear viscosity to bulk viscosity, and the ratio $Bq2=\kappa s/R\eta 1$ between surface dilation viscosity to bulk viscosity. Equation (6) becomes

introducing the surface Péclet number $Pes=UR/Ds$, a ratio of diffusion time, $R2/Ds$, to the characteristic time for particle motion, *R*/*U*. The force expression (17) remains unchanged, and the integral over *S _{p}* is split into parts that arise from the hemispheres in the two fluids with their corresponding $\pi 1$ and $\pi 2$. The constitutive Eq. (1) as well as the mass jump balance (9) and boundary condition (13) remains unchanged in their dimensionless form, while the initial and boundary conditions (7) and (8) for the surfactant concentration and velocity field become $\Gamma (r,t=0)=1$ for $r\u2208SI,\u2009\Gamma (r,t)=1$ for $r\u2208\u2202Sb,\u2009n\xb7[\u2207s\Gamma (r,t)]=0$ for $r\u2208\u2202Sp,\u2009u=ex$ for $r\u2208Sp$, and $u=0$ for $r\u2208Sb$.

## III. NUMERICAL METHOD

The problem is composed of two sub-problems: the Stokes equation in the fluid domain and the SCD equation at the interface. These two sub-problems are coupled and must be solved together. We employed the finite element method (FEM) in a domain of finite size. To avoid tracking the moving particle, the equations are solved in a frame that is moving with the particle, i.e., the particle is kept stationary at $x=(x,y,z)=0$, and the walls move with velocity $\u2212Uex$. Note that all our results are presented in a frame where the particle is moving and the walls are stationary. The Stokes equation (21) is thus solved for with a vanishing velocity on the particle surface $\u2202Sp$, the dimensionless boundary condition $u=\u2212ex$ on the outer surfaces of the simulation box, together with the stress boundary condition (22) at the fluid–fluid interface, using constitutive equation (23). The solution to the Stokes problem yields the velocity and pressure fields $ui(x)$ and $pi(x)$, hence $\tau i(x)$. To include the effect of the surfactant, in this work two approaches are employed: an implicit approach and an explicit approach. In the implicit approach, the SCD is solved together with the Stokes equation in one system. The non-linearity in the SCD is handled by a Picard iteration. In the explicit approach, the velocity at the fluid–fluid interface is used to solve the SCD Eq. (24) with boundary and initial conditions (7). After solving this equation, the new Γ is used to update interface tension *γ* using the interfacial equation of state, to be specified in Eq. (25). The explicit approach typically works well for problems at low Ma, i.e., for low interfacial elasticity. However, at larger values of Ma, the time step must be chosen extremely small for stable simulations. For those cases, the implicit approach is used, which is more expensive per time step, but which allows for an arbitrarily large time step. Steady-state results are obtained by running the simulation until the dimensionless time *t* = 100, which was found to be sufficiently long (see supplementary material Sec. S4).

We used an in-house FEM code to solve the full problem, i.e., the complete set of equations. A validation and convergence study is provided in the supplementary material Secs. S1(ii) and S1(iii). We make use of the open-source mesh generator Gmsh,^{44} which offers great control over the local element size in the domain (Fig. 2). Due to symmetry in the *xy*-plane, only half of the domain is meshed, and appropriate symmetry conditions are employed. In order to diminish boundary effects, a large simulation box is used: *L* = 3200 if not otherwise mentioned, which was shown to be sufficiently large (see supplementary material Sec. S2). The spherical particle is symmetrically immersed at the interface and translates with a constant velocity in the *x*-direction. The interface is located at *y* = 0, and it divides the simulation box into two domains: domain $\Omega 1$ in *y* < 0 is filled by a viscous liquid and $\Omega 2$ in *y* > 0 is filled by an inviscid fluid. The Stokes equation is solved in $\Omega 1$. Details of the numerical implementation are discussed by Dietrich *et al.*^{37} and Carrozza *et al.*^{45} We tested mesh convergence and used our mesh M3 as a compromise between quality and efficiency if not mentioned otherwise [supplementary material Sec. S1(iii)]. To solve Eq. (24), second-order finite difference schemes are used to approximate the time derivatives.

## IV. RESULTS AND DISCUSSION

In all the results presented in this work, we assume a linear interfacial equation of state: $\gamma =\gamma 0\u2212\Gamma kBT$ (dimensional form) where *γ*_{0} is the surface tension of a clean interface, and thus $\Pi s=K\pi =\Gamma kBT$ (dimensional form). Making the equation of state dimensionless using Eq. (20), we obtain

giving rise to the dimensionless capillary number $Ca=\eta 1U/\gamma 0$ (ratio between bulk stress and interface tension force) and Marangoni number $Ma=\Gamma 0kBT/\eta 1U$,^{46} the ratio between the force tending to deform the interface and surface elasticity, which tends to restore the original shape of the film. Note that from this we obtain $\Pi s=K\pi =Ma\u2009\Gamma $ in our dimensionless form, which shows that the Marangoni contribution to the force is proportional to Ma. Assuming a linear interfacial equation of state is a strong assumption, employed in most of previous works, and to be discussed later below.

The presented equations and methods capture the case of curved fluid-fluid interfaces. However, for this study, phase 2 is considered inviscid, therefore $\eta 2=0,\u2009\tau 2=0$, and the dimensionless number *λ* vanishes. Moreover, the interface is assumed to remain straight, i.e., the capillary number $Ca\u226a1$ and the mean curvature vanishes, which is a valid assumption for the typical particle sizes used in experiments.^{37} As the last term in Eq. (22) vanishes for a flat interface, results do not depend on the choice of Ca. Due to this assumed linear relationship between *γ* and Γ, we are left with dimensionless parameters $Bq1,\u2009Bq2$, Ma, and $Pes$, which completely describe the problem.

We begin by investigating the situation in the absence of any extra surface stresses, $Bq1=Bq2=0$, in Sec. IV A. Since one of our main objectives is to distinguish between the Marangoni and dilatational ($Bq2$) effects on the interface compressibility, we are going to introduce them separately here in Sec. IV B and to then study the interface when both mechanisms are present (Sec. IV C). It is worth noting at the outset that the case of vanishing $Bq2=0$ is an artificial case that has no analogue in the world of real interfaces. It is however important to delineate the contributions from shear and dilation to the measurable quantities and to allow for a comparison with previous results by Elfring *et al.*^{12} The ratio

which describes the viscous resistance to dilation, has been used^{9,11} instead of $Bq2$ to present and discuss results, as this ratio controls which deformation mode is “cheaper,” i.e., shear or dilation/compression. An analogy can be made to linear elastic materials, where the Poisson ratio can be expressed in terms of the shear and bulk (compression) modulus. Whenever we fix $Bq1$ and vary $Bq2$ or vice versa, our graphs describe the effect of Θ as well. The interface will be considered to exhibit a non-vanishing shear viscosity in Sec. IV D, while in Sec. IV E we focus on the influences of both interface viscosities on the drag coefficient. To guide a reader, the structure of the manuscript is further visualized in Table I.

Section . | $Bq1$ . | $Bq2$ . | Ma . | $Pes$ . | f
. | $us$ . | Figures . | Supplementary material figures . |
---|---|---|---|---|---|---|---|---|

III | ⋯ | ⋯ | ⋯ | ⋯ | ⋯ | ✓ | ⋯ | S2–S4 |

⋯ | ⋯ | ⋯ | ⋯ | ✓ | ⋯ | ⋯ | S5(a) | |

⋯ | ✓ | ⋯ | ⋯ | ✓ | ⋯ | ⋯ | S5(b), S6 | |

⋯ | ⋯ | ✓ | ✓ | ⋯ | ⋯ | ⋯ | S9 | |

IV A | ⋯ | ⋯ | ⋯ | ⋯ | ✓ | ✓ | ⋯ | S7,S8 |

⋯ | ⋯ | ✓ | ✓ | ✓ | ✓ | 3–,6 | S11,S12 | |

IV B | ⋯ | ✓ | ⋯ | ⋯ | ✓ | ✓ | 7–,8 | ⋯ |

IV C | ⋯ | ✓ | ✓ | ✓ | ✓ | ✓ | 9–,11(b) | ⋯ |

IV D | ✓ | ⋯ | ⋯ | ⋯ | ⋯ | ✓ | 12–,14 | ⋯ |

✓ | ⋯ | ✓ | ✓ | ⋯ | ✓ | 12 | ⋯ | |

✓ | ✓ | ⋯ | ⋯ | ⋯ | ✓ | 15 | ⋯ | |

✓ | ✓ | ✓ | ✓ | ⋯ | ✓ | 16 | ⋯ | |

✓ | ✓ | ✓ | ✓ | ⋯ | ⋯ | 17 | ⋯ | |

✓ | ✓ | ✓ | ✓ | ⋯ | ✓ | 18 | S13 | |

IV E | ✓ | ✓ | ✓ | ✓ | ✓ | ⋯ | 19 | ⋯ |

✓ | ✓ | ⋯ | ⋯ | ✓ | ⋯ | 20 | ⋯ | |

✓ | ✓ | (✓) | (✓) | ✓ | ⋯ | 21 | S10 |

Section . | $Bq1$ . | $Bq2$ . | Ma . | $Pes$ . | f
. | $us$ . | Figures . | Supplementary material figures . |
---|---|---|---|---|---|---|---|---|

III | ⋯ | ⋯ | ⋯ | ⋯ | ⋯ | ✓ | ⋯ | S2–S4 |

⋯ | ⋯ | ⋯ | ⋯ | ✓ | ⋯ | ⋯ | S5(a) | |

⋯ | ✓ | ⋯ | ⋯ | ✓ | ⋯ | ⋯ | S5(b), S6 | |

⋯ | ⋯ | ✓ | ✓ | ⋯ | ⋯ | ⋯ | S9 | |

IV A | ⋯ | ⋯ | ⋯ | ⋯ | ✓ | ✓ | ⋯ | S7,S8 |

⋯ | ⋯ | ✓ | ✓ | ✓ | ✓ | 3–,6 | S11,S12 | |

IV B | ⋯ | ✓ | ⋯ | ⋯ | ✓ | ✓ | 7–,8 | ⋯ |

IV C | ⋯ | ✓ | ✓ | ✓ | ✓ | ✓ | 9–,11(b) | ⋯ |

IV D | ✓ | ⋯ | ⋯ | ⋯ | ⋯ | ✓ | 12–,14 | ⋯ |

✓ | ⋯ | ✓ | ✓ | ⋯ | ✓ | 12 | ⋯ | |

✓ | ✓ | ⋯ | ⋯ | ⋯ | ✓ | 15 | ⋯ | |

✓ | ✓ | ✓ | ✓ | ⋯ | ✓ | 16 | ⋯ | |

✓ | ✓ | ✓ | ✓ | ⋯ | ⋯ | 17 | ⋯ | |

✓ | ✓ | ✓ | ✓ | ⋯ | ✓ | 18 | S13 | |

IV E | ✓ | ✓ | ✓ | ✓ | ✓ | ⋯ | 19 | ⋯ |

✓ | ✓ | ⋯ | ⋯ | ✓ | ⋯ | 20 | ⋯ | |

✓ | ✓ | (✓) | (✓) | ✓ | ⋯ | 21 | S10 |

### A. Surfactant dominant regime

While the drag coefficient on the particle had been studied before by Avramov *et al.*^{47} in the absence of extra surface stresses, we here present the most relevant results (including interface compressibility aspects) for this special case, as they serve as a benchmark for the results in Secs. IV B–IV E, and to validate our implementation.

When a particle translates at the interface covered with a surfactant, it pushes the surfactant, resulting in an accumulation of surfactant in the front of the particle, and a depletion at its rear. The surfactant concentration around a translating particle, for $Pes=1$, is shown in Fig. 3(a). This difference in the surfactant concentration causes a Marangoni flow; a flow from high to low concentration.^{48} Surfactant concentration profiles in the *x*-direction along a line at *z* = 0, passing through the center of the spherical particle at different Ma are shown in Fig. 3(b). To quantify the surfactant variation, we calculate the difference in the maximum and minimum concentrations, $\Delta \Gamma $. The inset in Fig. 3(b) shows $\Delta \Gamma $ at different Ma indicating that the relative surfactant variation is high for low Ma and decreases with increasing Ma.

The Marangoni flow is also highly affected by the $Pes$ (supplementary material Fig. S11). At very small $Pes\u226a1$, the time required for recovering homogeneity of the surfactant concentration due to fast diffusion of the surfactants relative to the particle motion is very short; therefore, uniformity is obtained instantaneously, $t\u226a1$. At $Pes<0.02$, there are negligible surfactant variations at the interface and at $Pes=0.001$ the surfactant concentration remains uniform. Beyond $Pes>1$, the variation does not change much upon further increasing $Pes$.

The variation in surfactant concentration, in turn, creates a gradient in $\u2207\gamma $ or $\u2207\Pi s$, which acts against the surface compression in front of, and dilation at the rear of the particle.^{49,50} Analogous to the bulk compressibility, the compressibility of monolayers can be defined as $C=\u2212\Gamma dA/d\Pi s|T$, where $A\u223c\Gamma \u22121$ is the area per surfactant molecule at *S _{I}*.

The interface symmetrically compresses at the front of the particle and dilates at the rear. The maximum dilation at the rear of the interface, denoted as $(\u2207s\xb7us)max$, serves as an indicator for the maximum interface compressibility. Barentin *et al.*^{16} argued that the monolayer can be considered incompressible if $|\u2207s\xb7us|\u226a1$, as one can write the incompressibility condition in the bulk close to a planar interface as $\u2207\xb7u=\u2207s\xb7us+\u2202uy/\u2202y=0$. Provided that $\u2202uy/\u2202y=0$, the compression $\u2207s\xb7u$ near the interface vanishes.

To further demonstrate the effect of coupling to the subphase on the compressibility of the clean interface, we varied the depth of the subphase. For a very shallow subphase there is a “layered” structure in the subphase (supplementary material Fig. S7), the normal velocity basically vanishes, and the interface becomes less compressible [supplementary material Fig. S8(a)]. Moreover, as the subphase depth is decreased, the bulk stresses, and thus the drag coefficient, diverge [supplementary material Fig. S8(b)]. These effects can be logically attributed to an increased importance of the subphase for decreasing subphase depth, inherently limiting the sensitivity of the probe to properties of the interface, as also concluded by Elfring *et al.*^{12} For a deep subphase, two different mechanisms can make an interface effectively incompressible: Marangoni effects and a viscous resistance to compressional flow.

Figure 4(a) shows the effects of $Pes$ on the interface compressibility at $Ma=10$. Maximum interface dilation $(\u2207s\xb7us)max$ for $Ma=10$ and $Ma=20$ as a function of $Pes$ are shown in Fig. 4(b). These results show that at the limits of small and high $Pes$, interface compression is independent of Ma. At $Pes\u2248O(1)$, an increase in Ma from 10 to 20 reduces the interface compressibility. These results show that when there is an increase in Ma in order to keep the interface compression constant, we should decrease $Pes$. Following the reasoning of Chisholm *et al.*, we can use the following simple argument to justify these results.^{51} Exploiting the interfacial equation of state, we can rewrite the SCD in the following form:

where the steady state is considered for simplicity. This equation shows that at high Ma the convection term is less important. The diffusion term at the right hand side also vanishes when $(Ma\u2009Pes)\u22121\u226a1$, which shows when there is an increase in Ma, in order to recover the interface compression $Pes$ should decrease. An increase in Ma also increases the Marangoni contribution to the force on the particle (in the regime where the linear equation of state holds), which consequently increases the drag coefficient, see Figs. 4(c) and 4(d). This is a known result first observed by Plateau,^{52} by comparing the damping of a magnetic needle immersed in a liquid with one placed at the liquid surface. He found that the resistance against deformation was higher at the surface than in the bulk and attributed this to a non-vanishing extra surface shear viscosity. Marangoni realized that trace contaminations are almost always present at liquid interfaces and that a non-uniform distribution of surface-active material, together with an interfacial equation of state, is responsible for Plateau's observation.^{53} Interestingly, 100 years after the works of Plateau and Marangoni the role of such impurities and Marangoni contribution in coalescence is still being resolved.^{54} It should be pointed out that such effects will be less strong when the surfactant concentration is higher, and the interfacial equation of state is no longer linear, but a weaker dependency of the surface tension with increasing concentration is typically observed.

To investigate the effects of Ma, we study the interface at different Ma at constant $Pes=1$. The interface compression and $(\u2207s\xb7us)max$ are shown in supplementary material Fig. S12. Results show that the interface compressibility decreases significantly by adding a small amount of surfactant to the interface. If $Ma>50$, its dependency on the concentration decreases, and there is a very slight change in the interface compression. As we have already pointed out, the Marangoni flow imposes a force on the particle and we expect the drag coefficient *f* to increase. Drag coefficients change from a diffusion dominant regime at $Pes=0.05$ to a convection dominant regime at $Pes=5$ [Fig. 5(a)]. The Marangoni contribution *F ^{M}* to the force amplitude as a function of Ma at different $Pes$ is shown in Fig. 5(b). The

*F*increases with Ma, and it is higher at high $Pes$. Moreover, Fig. 5(a) shows how the value of the drag coefficient converges to 11.5 for large $Pes$ and Ma. This value is consistent with that for an incompressible interface, as will be shown in Sec. IV B.

^{M}To explicitly show that surface dilation depends on both $Pes$ and Ma, we show the maximum interface dilation, normalized by the value for a clean, inviscid interface, as a function of both $Pes$ and Ma in Fig. 6. In this plot, the values on the axes are similar to the results for a clean, inviscid interface [supplementary material Sec. S1(ii)], since at small Ma the effect of Marangoni stresses becomes negligible, and at small $Pes$ diffusion dominates over the convection in the surface convection-diffusion (SCD) equation. Hence, in the latter case, the diffusion timescale for recovering a uniform distribution of the surfactant at the interface is much shorter than the convection timescale. In the limit of $Pes\u21920$, the surfactant concentration remains uniform, and regaining the uniform distribution is instantaneous. For values of Ma and $Pes$ that are both non-zero, the contourlines in this plot can be described well by a function of the form $Ma\u2009Pes=constant$, clearly showing that is the product between $Pes$ and Ma that plays a crucial role in surface dilation.

### B. Purely dilatational interfaces

For a purely expanding or compressing interface or for an interface with $Bq2\u226bBq1$, the interface stress tensor (2) simplifies to

Using the above relation in the tangential momentum jump at the interface, Eq. (22), we have

revealing the dilatational properties of the surface to be the surface tension *γ* and the surface dilatational viscosity $Bq2$. The expression in the parenthesis on the right hand side is also known as dynamic surface tension and can be expressed as

and is thus composed of both thermodynamic and viscous parts. Identifying the reversible and irreversible contributions to the dynamic surface tension is a difficult task, although it can be achieved by realizing that the thermodynamic part will be influenced by mass transport considerations (and hence geometry) and the viscous part only depends on the local strain rate.

For the limiting case of $Pes\u21920$, the surfactant remains uniform, the interface tension is constant, and the related term in Eq. (29) does not play any role for the dynamics at the interface. If we assume that the left hand side is a finite constant value, when *κ ^{s}* is getting larger, the interface compression $\u2207s\xb7us$ is getting smaller. Figure 7 shows our simulation results for the interface compression in the presence of the spherical particle that moves at constant speed in

*x*-direction. The maximum dilation $(\u2207s\xb7us)max$ is shown in the inset. We find that increasing $Bq2$ helps to make the interface more incompressible, while for $Bq2>1000$ the interface dilation does not change anymore significantly. The drag coefficient increases (by about one third) with increasing $Bq2$ (Fig. 8) and is expected to reach a constant value at $Bq2\u2192\u221e$, when the interface is in an incompressible state, as shown in Fig. 8,

*f*levels off at 11.5. As can be seen in Fig. 7, at this stage the interface is incompressible. These findings are consistent with the results for large $Pes$ and Ma, as shown in the previous section, but the origin of the surface incompressibility is different in both cases: here, rheological surface stresses cause a resistance against surface compression, whereas in Sec. IV A surface incompressiblity arises due to Marangoni stresses, which are thermodynamic in nature.

### C. Marangoni effects at dilatational interfaces

When $Pes\u2248O(1)$ or higher, both Marangoni flow and dilational viscosity can be responsible for the interface incompressibility. As we have shown in Eq. (29), it is not easy to distinguish these two phenomena. The relative importance of these two depends on $Bq2$ and Ma. Figure 9 shows the surfactant concentration profiles and $\Delta \Gamma $ at different $Bq2$. At high $Bq2$, the Marangoni flow diminishes due to the suppression of the surfactant gradient. These results indicate that the Marangoni flow is important at low $Bq2$ when the interface is still compressible. Maximum dilations $(\u2207s\xb7us)max$ and $\Delta \Gamma $ at different $Bq2=0,1,10,100$ as a function of surfactant concentration are shown in Fig. 10(a). At low surfactant concentration, the effect of $Bq2$ in both quantities is significant. Calculating the magnitude of the Marangoni contribution to the force also shows that *F ^{M}* decreases with increasing $Bq2$ [Fig. 11(a)].

At high $Bq2=100$, interface dilation only changes marginally upon addition of the surfactant to the interface. The result [Fig. 11(a)] indicates that the Marangoni effects are still present. We can consider the interface with $Bq2=100$ as being “almost incompressible,” and adding the surfactant makes the interface incompressible.

Comparison [Fig. 11(b)] of the drag coefficient at the interface with $Bq2=0$ and $Bq2=1$ shows negligible changes at different Ma. The sensitivity of *f* on Ma diminishes at high $Bq2$, where the interface is incompressible, and Marangoni effects absent [Fig. 11(b)].

### D. Interplay between surface viscosities, Marangoni flow, and (in)compressibility of the interface

So far, we have assumed that the interface has a vanishing shear viscosity. However, the presence of an interfacial shear viscosity will also lead to an energy dissipation when shear deformations occur. The detailed nature of the interface flow field hence depends on the relative contributions of other parameters such as interface compressibility and Marangoni flow at the interface. Yet, in microrheological experiments often only the shear viscosity is calculated from the observed drag forces or diffusivities and assumptions about incompressibility are typically made. It is thus worthwhile to investigate the interplay between these effects and the interfacial shear viscosity.

A compressible interface exhibits the possibility for the occurrence of Marangoni flow at the limit of $Bq2=0$. For an interface with $Pes\u21920$, the result shows that at $Bq1<10$, interface dilation increases with increasing $Bq1$ which means shear deformations cost more and more energy, the interface appears even more compressible, see Fig. 12. At $Bq1>10$, a further increase in surface shear viscosity slightly decreases the interface compressibility, most likely as, due to even a small shear deformations at these higher viscosities, the dissipation of momentum now is more localized close to the particle. In general, the interface shear viscosity does not have strong effects on the interface compressibility, as the effects are mainly indirect.

For a monolayer with a finite Péclet number $Pes=0.5$ (and $Ma=10$), for which now Marangoni effects are present, the interface dilation is also shown in Fig. 12. The behavior is similar to the case of $Pes\u21920$, where the surfactant remains uniform. Results show that with increasing the interface shear viscosity, the interface compressibility increases.

An illustration of an increase in the interface compressibility with $Bq1$ is suppression of the backflow in a confined flow. When a particle moves at the compressible interface, due to the confinement, it generates a backflow, where a part of the fluid in front of the particle moves in the opposite direction to compensate the motion of the particle. Figure 13 shows the interface velocity fields for low $Bq1=0.1$ and high $Bq1=100$ at $Pes\u21920$ and $Bq2=0$. The simulation domain length is *L* = 100. At low interface shear viscosity, there is a vortex at the interface. We can see that increasing the interface shear viscosity suppresses the vortex and causes the interface to move more uniformly. The center of the vortex, located at $(0,0,Zb)$, can be an indicator of the interface compressibility. With decreasing compressibility, the backflow increases and *Z _{b}* is getting closer to the particle surface and vice versa.

Profiles of the *x*-component of the surface velocity along the *z*-direction, from the particle surface, at different $Bq1$ are shown in Fig. 14, with the corresponding *Z _{b}* vs $Bq1$ in the inset. These results demonstrate that upon increasing the interface shear viscosity,

*Z*moves further away from the particle, which indicates an increase in the interface compressibility with increasing

_{b}*z*.

The uniform motion of the interface at high $Bq1$ can also be observed in an unbounded flow when walls are at infinity. Figure 15 shows velocity field, interface dilatation, and shear component of surface stress tensor for an interface with finite compressbility, $Bq2=0.1$, at three different interface shear viscosity $Bq1\u2208{0.1,1,5}$ and limit of $Pes\u21920$. The shear component of the surface velocity gradient is necessary to fully characterize the velocity gradient. The velocity fields show clearly the uniform motion at high $Bq1$. This uniform motion minimizes surface shear and increases the cost for dissipation of energy. This effect can be seen from the third column in Fig. 15. At low $Bq1$, the shear effect is very important at close distance from the surface of the particle, where there is a large gradient in shear component of the stress tensor. At high shear viscosity, the gradient decreases and shear effects happen at far distance from the particle.

Suppression of the vortex and uniform motion also affects the Marangoni flow. Figure 16 shows the velocity field $us$ and the surfactant concentration at $Ma=10$ and finite $Pes=0.1$. The interface has the same shear and dilatational viscosity as Fig. 15. To highlight the qualitative role of the Marangoni flow, its contribution to the surface velocity field, $\Delta us=us\u2212us0$ is also shown in Fig. 16, where $us0$ denotes the velocity field of the interface with uniform surfactant, at $Pes\u21920$. At low $Bq1$, the Marangoni flow occurs near the particle surface. Suppression of the vortex at higher shear viscosity enforces the surfactant to find another path to flow. Therefore, the flow is less focused and occurs even at far distance from the particle (see the second column). The results also show that the maximum flow is getting further from the particle with increasing $Bq1$. This is likely due to moving the vortex by increasing the interface compressibility. This difficulty for the surfactant to flow creates a higher gradient in front and back of the particle (third column), which consequently increases the Marangoni force acting on the particle.

The Marangoni contribution *F ^{M}* to the force and the interface viscous drag

*F*amplitudes are compared in Fig. 17. At low $Bq1$,

^{s}*F*increases significantly with $Bq1$, and it reaches a constant value with increasing $Bq1$. At high $Bq1$, the

^{M}*F*does not make a significant contribution to the total drag on the particle. An increase in

^{M}*F*with increasing surface shear viscosity was also reported by Elfring

^{M}*et al.*

^{12}for the motion of a disk at the interface. This can again be attributed to the enhancement of dilatation/compressional deformations as shear deformations become more costly.

When considering shear and dilational viscosities, a wide range of values has been reported (Table II), but for simple surface active components, such as small molecular weight surfactants, the ratio of shear and dilatational surface viscosities has been suggested to be of comparable magnitude, $\Theta =Bq2/Bq1=O(1)$.^{16,18} Depending on the chemical characteristics of the system at hand and the mobility of the surfactants (small molecular weight vs polymeric), the Péclet number might either be small or large compared to unity, we are going to discuss both types upon choosing $Pes\u2208{0.1,2.0}$. All qualitative effects are captured by a finite $Bq1=1$ and moderate $Ma=10$, which is the choice we are making for the present purpose.

Surfactant . | η_{1}
. | R
. | $\Gamma 0$ . | U
. | D
. _{s} | η
. ^{s} | κ
. ^{s} | $Bq1$ . | $Bq2$ . | Ma . | $Pes$ . |
---|---|---|---|---|---|---|---|---|---|---|---|

(m Pa s) . | (μm)
. | (10^{6}/μm^{2})
. | (μm/s)
. | (μm^{2}/s)
. | (Ns/m) . | (Ns/m) . | . | . | . | . | |

DPPC | 0.89 | 170 | 0.98 | 500 | 70 | $10\u22125$ | 1 | 66 | 6.6 $\xd7106$ | 9200 | 1200 |

F-108 | 0.89 | 1 | 0.17 | 1 | 30 | 0 | 1.36 | 0 | 1.53 $\xd7109$ | 8 $\xd7105$ | 0.04 |

P-123 | 0.89 | 1 | 0.25 | 1 | 20 | 0 | 0.17 | 0 | 1.96 $\xd7108$ | 1.2 $\xd7106$ | 0.03 |

H_{21}T_{8}H_{21} | 0.89 | 1 | 0.32 | 1 | 170 | $2.5\xd710\u221211$ | $2.6\xd710\u221211$ | 38 | 83.7 | $1.2\xd710\u22122$ | 250 |

Surfactant . | η_{1}
. | R
. | $\Gamma 0$ . | U
. | D
. _{s} | η
. ^{s} | κ
. ^{s} | $Bq1$ . | $Bq2$ . | Ma . | $Pes$ . |
---|---|---|---|---|---|---|---|---|---|---|---|

(m Pa s) . | (μm)
. | (10^{6}/μm^{2})
. | (μm/s)
. | (μm^{2}/s)
. | (Ns/m) . | (Ns/m) . | . | . | . | . | |

DPPC | 0.89 | 170 | 0.98 | 500 | 70 | $10\u22125$ | 1 | 66 | 6.6 $\xd7106$ | 9200 | 1200 |

F-108 | 0.89 | 1 | 0.17 | 1 | 30 | 0 | 1.36 | 0 | 1.53 $\xd7109$ | 8 $\xd7105$ | 0.04 |

P-123 | 0.89 | 1 | 0.25 | 1 | 20 | 0 | 0.17 | 0 | 1.96 $\xd7108$ | 1.2 $\xd7106$ | 0.03 |

H_{21}T_{8}H_{21} | 0.89 | 1 | 0.32 | 1 | 170 | $2.5\xd710\u221211$ | $2.6\xd710\u221211$ | 38 | 83.7 | $1.2\xd710\u22122$ | 250 |

The surface flow field and its characteristics, as well as the surfactant concentration profiles, are shown for three different $\Theta \u2208{0.1,1,10}$ and $Pes=0.1$ in Fig. 18 (corresponding data for $Pes=2.0$ are available in supplementary material Fig. S13). The different components of the velocity field in Fig. 18 again reveal the interplay between the shear and dilational/compression. This figure demonstrates how the relative contributions of compression/dilation and the shear deformations depend on the relative values of Θ and the magnitudes of the Marangoni and bulk to interfacial stresses, congruent with work on macroscopic rheometers, which employ mixed flows, such as the double wall sinusoidal ring^{10} and the interfacial Cambridge extensional rheometer.^{11,57}

### E. Influences of the surface viscosities on the drag coefficient of a spherical particle

Although recent available techniques make the experimental realization of surface microreology relatively straightforward, especially in the particle tracking techniques, one has to rely on a hydrodynamic model of the monolayer to obtain interface quantities such as the interface shear viscosity. Particle tracking techniques are a branch of, so-called, microrheological techniques, which have advantages over macrorheology including smaller sample sizes, access to higher frequency ranges, the ability to measure small-scale material heterogeneities, and greater sensitivity.^{9} In particle tracking, the position of a Brownian probe is recorded as a function of time, and the diffusivity *D* of the probe is calculated from a mean squared displacement (MSD) which in two dimensions is related to diffusivity by

where *r* is the distance between a pair of particles, *τ* is the lag time, and and $\u3008\u3009$ denotes an ensemble average. In general, *D* is related to the particle's drag coefficient *f* by the Einstein-Sutherland relation^{58,59}

where $kB$ is Boltzmann's constant and *T* is the absolute temperature. To calculate the surface shear viscosity using (32), it is necessary to investigate the proposed model for the drag coefficient as a function of the available monolayer parameters.

The drag coefficient vs $\Theta =Bq2/Bq1$ for both $Pes$ (in Figs. 18 and S13) is shown in Fig. 19. For the larger $Pes$, both the drag coefficient and the flow fields are basically unaffected by Θ over the experimentally relevant range. The surface flow field is only weakly compressed and dilated, while the surfactant concentration is distorted by about 10%. At the smaller $Pes$ the situation is reversed (Fig. 18). While the surfactants remain homogeneously spread at the interface, the surface velocity field exhibits larger shear gradients and weaker compression/dilation effects with increasing Θ.

According to these results, the effect of $Pes$ may be best detected from the Γ variation, while Θ affects the *z*-component of the surface velocity field, evaluated along the *x* =* z*–diagonal. The above discussion pertains to small molecule surfactants;^{12,60} for other types of structure interfaces, the Boussinesq–Scriven model is probably not adequate and strong viscoelastic effects are present, and dilatational elasticity becomes important.^{9} However overall the trends and ideas observed in the simulations can be used to assess under which conditions the drag force is going to be determined by shear viscosities.

While Fig. 19 represents the drag coefficient for an interface in general cases when all dimentionless parameters are present, models for the drag coefficient of protruding particle as a function of interface shear viscosity have been reported for the limiting cases of compressible and incompressible interfaces. Danov *et al.*^{18} solved numerically the problem of a translating particle at a compressible Newtonian interface with constant surface tension. They reported the particle diffusion coefficient $\u221d1/f$ at different contact angles as a function of $Bq1$ and $Bq2$. We simulated the interface under similar conditions. The drag coefficient *f* and magnitude of bulk *F ^{b}* and interface viscous

*F*drags are shown in Fig. 20. Comparing our results with Danov's study for a particle with contact angle $90\xb0$ shows that these two results are in a very good agreement [Fig. 20(b)], given that the numerical methods differ and that we are operating in a differently bounded domain. Bonales

^{s}*et al.*

^{61}have used Danov's theory to calculate the shear viscosity of two polymer Langmuir films. They compared the results with the data obtained by canal viscosimetry. Maestro

*et al.*

^{4}used this theory at the glass transition and compared the results with particle tracking techniques. Sickert

*et al.*

^{22}used Danov's theory in the study of monolayers of fatty acids and phospholipids. They have shown that the surface shear viscosity calculated from Danov's theory was lower than experimental data. A critical assessment of these findings can be found in Samaniuk and Vermant.

^{9}

The incompressibility of interfaces is a widely accepted assumption and has been used in many pioneering works in the problem of particle motion at the interfaces and membranes in both non-protruding particles^{3,13,15} and protruding particle.^{23,62} Fischer *et al.*^{23} solved the problem of translation of a particle with radius *R* in an incompressible viscous monolayer with the surface shear viscosity *η ^{s}*, between two viscous phases with viscosities

*η*

_{1}and

*η*

_{2}. They solved the incompressible surface Stokes equations

where $Fs$ is an external surface force parallel to the monolayer and $C|s$ represents the jump in *C* across interface. They present solutions to their model for drag coefficients on spherical particles at interfaces for the limiting cases where $Bq1\u226a1$ or $Bq1\u226b1$. They have obtained the following result for the translational drag coefficient *f* as a series expansion of Boussinesq number $Bq=\eta s/(\eta 1+\eta 2)R$ for $0<Bq\u226a1$. For our setup $Bq=Bq1$ and thus

Fitted expressions for *f* from numerical results give the following formulas for *f*_{0} and *f*_{1}:^{23}

where *d* is the signed distance from the apex of the sphere to the plane of the interface, which defines the contact angle. Equation (35) had been published by Fischer *et al.* with $\alpha 0=32$. This value is however incompatible with the values for *f*_{0} mentioned in both their text and data. Upon adjusting *α*_{0} in (35) to exactly recover their published fit function, we obtain the corrected value $\alpha 0=36.3$ that should be used when making use of Eq. (35) (details in supplementary material Sec. S5, a more stringent test of the whole fitting formula would require studying ellipsoidal or spheroidal particles instead of spheres). Notice that when $d\u2192\u221e$, the expression inside the radical goes to 1 and $f0=6\pi $ and $f1=0$, which is the correct theoretical value for a sphere in the bulk (Stokes law). The two expressions for *f*_{1} do not converge for $d\u21920$ since the series expansion in Bq of Oseen's tensor no longer converges for $Bq\u21920$.^{23} For a translational drag on a half immersed sphere in a non viscous monolayer, they find $f\u2248f0\u224811.7$, which is higher than at a free surface ($f=3\pi $). This is due to surface incompressibility, which they attribute to Marangoni effects. The value of $f0\u224811.7$ obtained by Fischer *et al.* (the above fit formula with unrevised $\alpha 0=32$ instead gives $f0=11.08$) for an incompressible and inviscid interface is slightly higher than the value of $f\u224811.5$ we found for an incompressible interface due to Marangoni flow in Sec. IV A or dilational viscosity in Sec. IV B or combinations of both effects in Sec. IV C.

These results are important and show that the drag coefficient of a particle at incompressible interfaces is independent of the origin of the incompressibility (dilatational viscosity, Marangoni effects or a combination of both) and that its higher value can not only be related to the Marangoni effects, as suggested by Fischer *et al.* The conditions required for surface incompressibility were not investigated in detail by Fischer *et al.*^{23} Figure 11(a) shows that at $Bq2=100$ the Marangoni effects are still present but they do not have significant effects on the drag coefficient [Fig. 11(b)]. These results do not support their claim that the Marangoni effects become very noticeable in the limit of vanishing surface compressibility.^{23}

We compared our results with Fischer's study for $Bq1\u226a1$ for an incompressible interface with $Bq2=106$ using a linear fit to the data. At $Bq1\u21920$, the drag coefficient $f=f0\u224811.52$ is smaller than the value reported by Fischer *et al.,*^{23} and for the slope *f*_{1} we obtain $f1\u22486.85$. Drag coefficients of the particle for intermediate $Bq1$ ($1\u2264Bq1\u226410$) and high values ($10<Bq1\u2264104$) at an incompressible interface (with $Bq2=106$) are shown in Fig. 21. Our data can be fitted with linear equations,

as shown in the figure. The drag coefficients for the intermediate $Bq1$ were not reported by Fischer *et al.*, since there is no solution for that regime. A comparison between the drag coefficient at small and intermediate $Bq1$ clearly shows a noticeable difference between these possibly experimentally relevant two regimes (Table II).

A similar analysis for a compressible interface ($Bq2=0.1,Pes\u21920$) is captured by the analytical expressions $f=11.234+\u20091.931Bq1$ and $ln\u2009f=1.07+0.857\u2009ln\u2009(Bq1)$ for the intermediate and high $Bq1$, respectively, in the same range of $Bq1$ used for the incompressible fluid in Fig. 21.

Fischer's model has been used for the regime of $Bq1\u226a1$ to calculate the interface shear viscosity from *D*, with *R*, *η*_{1} and *T* at hand.^{4,9,22,30,63,64} However, if one utilizes it within the $0.1<Bq1<10$ regime, there is at least an order-of-magnitude agreement with estimates of *η ^{s}* calculated from an extrapolation of the two limiting solutions of the Fischer model.

^{9}The model presented in this work for drag coefficient at $1\u2264Bq1\u226410$—in dimensional form $f=\eta 1R(14.830+3.677\u2009Bq1)$—can be used for such an intermediate regime.

## V. CONCLUSIONS

This work offers a comprehensive study of the effects of interfacial properties on the drag of a spherical probe partially immersed into a viscous bulk fluid. We calculated the force on a spherical probe particle translating at constant velocity and centered at the flat interface, while the interface was initially covered uniformly by a surface active component. The translating particle compresses the interface at its front and dilates at its rear. We demonstrated how the quantitative occurrence of this phenomenon is affected by the properties of the surface active component, namely the concentration, the interfacial equation of state and the Péclet number, and the interface shear and dilatational viscosities. The presence of a surface active material was shown to reduce the interface compressibility and this effect strongly depends on the surface Péclet number. We quantified how a non-vanishing interface dilatational viscosity helps to reduce the local interface compressibility, while shear viscosity has only a smaller, opposite, and obviously indirect effect. A low value for the shear viscosity is found to provide an easy way to deform locally and reduces the divergence of the surface velocity field.

For the special case of inviscid interfaces carrying surfactants, the product between Ma and $Pes$ plays a crucial role in determining the dilation of the interface, and we show that for $MaPes>15$, the interface dilation is only 20% as compared to a clean, inviscid interface. Increasing the product between $Ma\u2009Pes$ further can thus make the interface effectively incompressible.

When interfacial viscosities are included as well, the Marangoni contribution *F ^{M}* to the forces experienced by the particle is found to be affected by the magnitude of both interfacial shear and dilatational viscosities. The particle induces a flow field that compresses/dilates the interface, which leads to surfactant concentration variations, and thus to variations in the interfacial tension along the interface. The flow field depends on the rheological properties of the interface, here to be taken purely viscous, and the interface Péclet number. A non-uniform distribution of surfactants leads to Marangoni flows, but these flows are coupled to the surface viscosities of the interface, resulting in a complex interplay of effects. Markedly, our simulation results showed that interface viscosities have an opposite effect on the Marangoni force

*F*acting on the particle; dilatational viscosity decreases but shear viscosity slightly increases

^{M}*F*.

^{M}The present work demonstrated that an interface can be incompressible by dilatational viscosity, Marangoni effects, or combinations of both. The drag coefficient of the particle at incompressible interfaces for these three cases, within the limit of vanishing interface shear viscosity, exhibited the same value of 11.5. We also calculated the particle's drag coefficient over a wide range of dimensionless interface shear viscosities ($Bq1$) for compressible and incompressible interfaces and provided approximate analytical expressions for them.

The results presented in this work validate conclusions from the literature that were limited to a narrow asymptotic regime and only valid for non-protruding particles in shallow subphases.^{12} The results of the present study should therefore help interpreting shear and dilatational rheological, as well as particle tracking experiments, where particles are generally protruding and the subphase should be deep in order to effectively probe the interface. Moreover, our numerical method allowed us to explore the parameter space far from the asymptotic regimes.

Here, we focused on a spherical probe symmetrically immersed at the interface, with a contact angle equal to $90\xb0$. An extension to non-spherical probes to study the effects of the particle aspect ratio, orientation, and contact angle is straightforward and subject of a subsequent study. We here assumed that the surface tension decreases linearly with increasing surfactant concentration. This assumption is limited to dilute condutions^{21,65} and should be released within the otherwise unchanged framework if one aims for a quantitative comparison with experimental data, provided the nonlinear features of the surface tension have been fully resolved beforehand, which would extend the linearized equation used here. Nonequilibrium molecular dynamics and Monte Carlo simulations along the lines indicated in recent works^{56,66–68} for complex interfaces may help to provide a route to calculate the missing ingredients, the material functions *γ,* and surface viscosities as a function of Ma for atomistically detailed or coarse-grained models of amphiphilic surfactants.

## SUPPLEMENTARY MATERIAL

See the supplementary material for the validation of the FEM program, mesh convergence, effects of subphase depths, transient behavior, the coefficient *α*_{0}, and additional data, table, and figures referred to within the manuscript.

## DEDICATION

This manuscript is dedicated to the memory of Dr. Robert “Bob” Byron Bird.

## ACKNOWLEDGMENTS

The authors thank M. A. Hulsen at the Eindhoven University of Technology (TU/e) for access to the TFEM software libraries, and A. Moghimikheirabadi for help with the compiling of Table II and Table S2. This project has been supported by the Swiss National Science Foundation through Grant Nos. 200021L_185052 (M.K.) and 20021_165974 (J.V.). M.K. acknowledges support from the Swiss National Supercomputing Centre (project s987).

## DATA AVAILABILITY

The data that support the findings of this study are available within this article.