A fine balance of chemotactic and hydrodynamic torques: when microswimmers orbit a pillar just once

TL;DR

This study investigates how self-propelling droplet microswimmers interact with microfluidic pillars, revealing that chemical and hydrodynamic forces cause them to orbit once before detaching, with implications for understanding microswimmer behavior.

Contribution

The paper introduces a combined experimental and modeling approach to quantify chemical and hydrodynamic torques influencing microswimmer detachment from pillars.

Findings

Droplets tend to orbit once before detaching due to chemical and hydrodynamic effects.

The Langevin dynamics model accurately describes swimmer behavior near pillars.

Chemical and hydrodynamic parameters depend on boundary conditions like wall curvature.

Abstract

We study the detention statistics of self-propelling droplet microswimmers attaching to microfluidic pillars. These droplets show negative autochemotaxis: they shed a persistent repulsive trail of spent fuel that biases them to detach from pillars of a certain size after orbiting them just once. We have designed a microfluidic assay recording swimmers in pillar arrays of varying diameter, derived detention statistics via digital image analysis and interpreted these statistics via the Langevin dynamics of an active Brownian particle (ABP) model. By comparing data from orbits with and without residual chemical field, we can independently estimate quantities like hydrodynamic and chemorepulsive torques, chemical coupling constants and diffusion coefficients, as well as their dependence on boundary conditions like wall curvature.

Equations

$$\dot{\mathbf{r}}$$

$$\dot{\mathbf{e}}$$

$$\dot{r}$$

$$\dot{\theta }$$

$$f(t_{\text{d}};\mu ,\lambda )=\sqrt{\frac{\lambda }{2\pi t_{\text{d}}^3}}\exp \left[\frac{-\lambda (t_{\text{d}}-\mu {)}^2}{2{\mu }^2{t}_{\text{d}}}\right],$$

$$f(t_d)=\sqrt{\frac{\lambda }{2\pi t_{\text{d}}^3}}\exp \left[\frac{-\lambda }{2t_{\text{d}}}\right].$$

$${\partial }_{t}c=D_{\text{m}}{\partial }_r\left(r{\partial }_{r}c\right)/r,$$

$$\text{MSD}=\{\begin{array}{l}4D_{\text{T}}\Delta t+v^2\Delta t^2\text{for}t\ll {\tau }_{\text{R}}\\ (4D_{\text{T}}+v^2{\tau }_R)\Delta t\text{for}t\gg {\tau }_R\end{array},$$

$$\delta Y=\sqrt{\sum _{i=0}^n{\left(\frac{\partial Y}{\partial x_i}\delta x_i\right)}^2}\text{ for }Y=Y(x_0,\dots ,x_n)$$

$$\frac{d\mathbf{e}}{dt}=\mathbf{\Xi }\times \mathbf{e},\text{with }\mathbf{\Xi }=\mathbf{\Omega }+\xi .$$

$${\partial }_{t}P(\mathbf{e},t)=-\sum _{i=1}^n{\partial }_{e_i}⟨\delta (\mathbf{e}-{\mathbf{e}}_0)\frac{d{e}_i}{dt}⟩.$$

$${\partial }_t⟨\delta (\mathbf{e}-{\mathbf{e}}_0)⟩$$

$$=-{\nabla }_e⟨\frac{d\mathbf{e}}{dt}\delta (\mathbf{e}-{\mathbf{e}}_0)⟩$$

$$=-{\nabla }_e⟨(\mathbf{\Xi }\times \mathbf{e})\delta (\mathbf{e}-{\mathbf{e}}_0)⟩$$

$$=-{\nabla }_e⟨(\mathbf{\Omega }\times \mathbf{e})\delta (\mathbf{e}-{\mathbf{e}}_0)⟩-{\nabla }_e⟨(\xi \times \mathbf{e})\delta (\mathbf{e}-{\mathbf{e}}_0)⟩$$

$$⟨\xi (t)R[\xi ]⟩=\underset{-\infty }{\overset{+\infty }{\int }}d{t}^{\prime }⟨\xi (t)\xi (t^{\prime })⟩=⟨\frac{\delta R[\xi ]}{\delta \xi (t^{\prime })}⟩,$$

$${\partial }_{t}P(\mathbf{e},t)=D_{\text{R}}{\text{L}}_e^{2}P(\mathbf{e},t)-{\nabla }_e\cdot (\mathbf{\Omega }\times \mathbf{e})P(\mathbf{e},t),$$

$${\partial }_{t}P(\mathbf{e},t)=O_{p}P(\mathbf{e},t).$$

$$P(\mathbf{e},t)=⟨\mathbf{e}|e^{D_{\text{R}}{L}^{2}t}|{\mathbf{e}}_0⟩.$$

$$P(\theta ,t|{\theta }_0,t_0)=\frac{1}{\sqrt{4\pi D_{\text{R}}t}}\exp (-\frac{(\theta -{\theta }_0{)}^2}{4D_{\text{R}}t}).$$

$$P(\theta ,t|{\theta }_0,t_0)=\frac{1}{\sqrt{4\pi D_{\text{R}}t}}\left[\exp (-\frac{(\theta -{\theta }_0{)}^2}{4D_{\text{R}}t})-\exp (-\frac{(\theta -(2{\theta }_{\text{e}}-{\theta }_0){)}^2}{4D_{\text{R}}t})\right].$$

$$S(t)$$

$$=\frac{1}{2}\left[\mathrm{erf}(\frac{{\theta }_{\text{e}}-{\theta }_0}{\sqrt{4D_{\text{R}}t}})+1-\mathrm{erf}(\frac{{\theta }_0-{\theta }_{\text{e}}}{\sqrt{4D_{\text{R}}t}})-1\right]$$

$$=\frac{1}{2}\left[\mathrm{erf}(\frac{{\theta }_{\text{e}}-{\theta }_0}{\sqrt{4D_{\text{R}}t}})+\mathrm{erf}(\frac{{\theta }_{\text{e}}-{\theta }_0}{\sqrt{4D_{\text{R}}t}})\right]$$

$$=\mathrm{erf}(\frac{{\theta }_{\text{e}}-{\theta }_0}{\sqrt{4D_{\text{R}}t}}).$$

$$f(t)$$

$$=\frac{({\theta }_{\text{e}}-{\theta }_0)}{\sqrt{4\pi D_{\text{R}}{t}^3}}\exp (-\frac{({\theta }_{\text{e}}-{\theta }_0{)}^2}{4D_{\text{R}}t}).$$

$$S(t)=$$

$$=$$

$$+\mathrm{erf}(\frac{{\theta }_0+{\theta }_{\text{e}}-2{\theta }_{\text{r}}}{\sqrt{4D_{\text{R}}t}})-\mathrm{erf}(\frac{{\theta }_0-2{\theta }_{\text{e}}+{\theta }_{\text{r}}}{\sqrt{4D_{\text{R}}t}})].$$

$$f(t)=$$

$$-({\theta }_0-{\theta }_{\text{e}})\exp (-\frac{({\theta }_0-{\theta }_{\text{e}}{)}^2}{4D_{\text{R}}t})+({\theta }_0-2{\theta }_{\text{e}}+{\theta }_{\text{r}})\exp (-\frac{({\theta }_0-2{\theta }_{\text{e}}+{\theta }_{\text{r}}{)}^2}{4D_{\text{R}}t})$$

$$+({\theta }_0+{\theta }_{\text{e}}-2{\theta }_{\text{r}})\exp (-\frac{({\theta }_0+{\theta }_{\text{e}}-2{\theta }_{\text{r}}{)}^2}{4D_{\text{R}}t})-({\theta }_0-{\theta }_{\text{r}})\exp (-\frac{({\theta }_0-{\theta }_{\text{r}}{)}^2}{4D_{\text{R}}t})].$$

$$f(t)=D_{\text{R}}\frac{({\theta }_{\text{e}}-{\theta }_0)}{\sqrt{4\pi {\bar{t}}^3}}\exp (-\frac{({\theta }_{\text{e}}-{\theta }_0{)}^2}{4\bar{t}}).$$

Full text

Biological microswimmers operate in complex geometries and react to chemical and physical gradients, both external, e. g. in nutrients (chemotaxis Madigan and Martinko (2005); Adler (1966); Hazelbauer (2012)), and self generated ones (autochemotaxis Bonner and Savage (1947); Zhao et al. (2013); Reid et al. (2012)), the latter enabling them to communicate and cooperate. The boundary conditions imposed e.g. by soil packings and interfaces affect both the active dynamics of swimmers, via hydrodynamic wall interactions Spagnolie and Lauga (2012); Berke et al. (2008) or active accumulation Li and Tang (2009); van Teeffelen and Löwen (2008); Elgeti and Gompper (2013), as well as the diffusive spread of chemotactic fields. The coupling of these dynamics leads to complex feedback phenomena that can significantly influence behaviour like arrest, aggregation and biofilm formation Tuval et al. (2005).

In order to derive and test tractable models for systems featuring such complex feedback it is important to develop experimental assays to decouple and control individual contributions.

We have previously demonstrated one quite unexpected phenomenon in a highly tunable artificial model system, while studying the interaction of self-propelling droplet swimmers with pillars in quasi-2D microfluidic cells. These droplets shed a chemorepellent that biases them towards detachment after circling the pillar once Jin et al. (2018). Our system has the advantage of directly comparing wall attachment with and without chemotactic repulsion, as well as controllable curvature set by the pillar size.

In this paper we present a quantitative study of the pillar interaction, analyzing the data from multiple experimental data sets in the context of an analytical Langevin model incorporating hydrodynamic, chemotactic and stochastic torques. From our quantitative analysis we gain insight into the respective strengths of wall attraction, chemorepulsion, as well as rotational diffusion and their effects on interfacial capture.

Our swimmers consist of oil droplets of radius $a\approx$50\text{,}\mathrm{\SIUnitSymbolMicro m}$$, dissolving gradually in a micellar aqueous surfactant solution, on a time scale of up to several hours (See SI Sec. I for experimental protocols). The oil solubilises by diffusing into surfactant micelles in a boundary layer around the droplet, a process that also depletes the surfactant coverage of the droplet’s oil-water interface Herminghaus et al. (2014); Maass et al. (2016). If the droplet position fluctuates, the symmetric cloud of filled micelles is deformed towards trailing behind the droplet, such that there are more empty micelles at the anterior. Empty micelles can disintegrate to replenish the surfactant coverage at the interface: thus, a moving droplet has more surfactant at its anterior, leading to a self-sustaining gradient in interfacial tension towards the posterior and a force propelling the droplet ($v\approx$ 30–60$\text{\,}\mathrm{\SIUnitSymbolMicro m}\text{\,}{\mathrm{s}}^{-1}$).

Chemotaxis and autochemotaxis follow naturally from the droplets’ mechanism of propulsion Jin et al. (2017). Droplets are attracted by empty micelles, such that they follow surfactant gradients (chemotaxis) and avoid areas of filled micelles, which act as a chemorepellent. The micelles diffuse slowly in comparison to the droplet motion, with a diffusion constant $D_{\text{m}}\approx\,$100\text{,}{\mathrm{\SIUnitSymbolMicro m}}^{2}\text{,}{\mathrm{s}}^{-1}$$ Candau et al. (1984), such that the repulsive trail of filled micelles persists over cruising ranges of tens to hundreds of droplet diameters. We have experimentally confirmed Jin et al. (2017) that this negative autochemotaxis is indeed mediated by the diffusion of filled micelles.

A persistent swimmer will not successfully attach to a convex wall without an attractive torque, simply for reasons of geometry. Since the required torque increases with curvature, swimmers experiencing a finite wall attraction are expected to scatter off small pillars and get trapped at large ones Sipos et al. (2015); Takagi et al. (2014); Simmchen et al. (2016). Our specific swimmers are attracted to walls and indeed exhibit scattering from small and trapping at large pillars. However, at pillars of intermediate radii, we observe a significant increase of detachment after circling the pillar once, i.e. when swimmers approach and re-enter their own trajectory Jin et al. (2018): the filled micelles diffuse away from the pillar (Fig. 4 (c)), providing an additional chemorepulsive gradient and favouring detachment. We illustrate this by consecutive video stills of a typical interaction with a $100\text{\,}\mathrm{\SIUnitSymbolMicro m}$ pillar in Fig. 2 and a histogram of detention lengths for many similar interactions in Fig. 3, which clearly shows a peak around lengths corresponding to one orbit (see also supporting movie, Sec. VII).

The long-tailed peaks in the statistics of detention times suggest stochastic influences on the pillar interaction Schaar et al. (2015); Spagnolie et al. (2015). Rotational noise and convex curvature cause a fraction of swimmers to detach almost immediately (Fig. 3, blue peak). For droplets attached long enough to approach their own residual chemical field, we can distinguish two cases: If the chemical repulsion is too strong to allow the droplets to re-enter their trajectory, this creates a non-crossing peak in the histogram of detention lengths slightly before one full orbit (green). Otherwise, the droplet is forced to swim along its own trail, as traced in Fig. 2. However, the added chemotactic repulsion adds a bias towards detachment, resulting in a self-crossing peak in the histogram of detention lengths slightly after one orbit (yellow). We note that we exclude data from non-crossing interactions from our analysis, since their chemical interactions break radial symmetry and would require rather sophisticated modelling. We can already estimate the forces and torques of interest from the following simple analytical model applied to the cases of short-time scattering and self-crossing.

For a quantitative estimate of the effects of curvature, wall interaction and autochemotaxis we model our swimmers as active Brownian particles (ABP) Romanczuk et al. (2012); Schmitt and Stark (2016). The swimmer is represented by a point particle moving with speed $v$ in the direction $\mathbf{e}$. The dynamics of position $\mathbf{r}$ and direction $\mathbf{e}$ are described by the overdamped Langevin equations

[TABLE]

$\bm{\Omega}$ is the sum of deterministic drifts of $\mathbf{e}$, and the rotational fluctuations are modeled by the rotational diffusivity $D_{\text{R}}$ and a normalized zero-mean Gaussian white noise $\bm{\xi}$ with variance $\langle\bm{\xi}_{i}(t)\bm{\xi}_{j}(t^{\prime})\rangle=\delta_{ij}\delta(t-t^{\prime})$.

We note that in active emulsions, the orientation can be affected by both whole body rotations and shifts in interfacial tension. However, in our reduction to ABP point particles, these revert to effective forces and torques.

The droplet motion is restricted to the $xy$ plane, with a perpendicular $\bm{\Omega}=(0,0,\Omega)^{T}$. In our system $\Omega=\Omega_{\text{W}}+\Omega_{\text{C}}$ comprises the sum of deterministic drifts of $\mathbf{e}$: the external (mostly hydrodynamic) torque due to wall attraction, $\Omega_{\text{W}}$, and autochemotaxis, $\Omega_{\text{C}}$, from the second orbit. A schematic representation of the system is shown in Fig. 4(a). We use polar coordinates $(r,\alpha)$ with the origin $r=0$ at the center of a pillar of radius $R$, with a swimmer of radius $a$ at the position $\mathbf{r}$ from the origin. We project the orientational unit vector of the swimmer $\mathbf{e}$ on two orthogonal unit vectors, a radial one pointing towards the pillar center, $\mathbf{e}_{\perp}$, and one tangential to the pillar interface, $\mathbf{e}_{\parallel}$. We define $\theta\in[-\pi/2,\pi/2]$ as the angle between $\mathbf{e}$ and $\mathbf{e}_{\parallel}$, with $\theta\approx 0$ when the droplet is moving along the pillar tangent. Without loss of generality, we assume the swimmer to move counterclockwise around the pillar, such that it approaches with $\theta>0$ and detaches with $\theta<0$. From Eq. (1), we derive Spagnolie et al. (2015)

[TABLE]

where $\omega>0$ is the angular velocity of the swimmer around the pillar. The swimmer detaches when $\theta$, starting from $\theta_{0}$, exceeds a fixed escape angle $\theta_{\text{e}}$ Zöttl and Stark (2016); Spagnolie et al. (2015) (Fig. 4(b)). The detention time $t_{\text{d}}$ is therefore the first-passage time of $\theta_{\text{e}}$ of a stochastic process described by Eq. (2). Here we take $\theta_{0}=0$, and $\theta_{\text{e}}=-0.962$ as the escape angle measured at the moment of detachment (see SI Sec. II.1).

Our restriction to counterclockwise motion would in principle necessitate a reflective boundary condition for $\theta_{r}=\pi/2$, however, for moderate rotational diffusion $D_{r}$ (see also SI Sec. II.6 and V), this head-on orientation is hydrodynamically unlikely in the case of our highly persistent pusher-type swimmers Herminghaus et al. (2014); Krüger et al. (2016) even in the limiting case of zero wall curvature and no chemotactic wall repulsion, and the boundary condition can be safely neglected.

Generally, $\Omega_{\text{W}}$ and $\Omega_{\text{C}}$ are functions of $r(t),\,\theta(t)$ and $\alpha(t)$, depending on the exact nature of the hydrodynamic and chemical fields. However, we can treat them as constants using the following approximations. We assume the transition between attached and detached state to be instantaneous compared to the detention time, such that $r\approx R+a$, $\theta\approx 0$ and $\Omega_{\text{W}}$ is constant for the swimmer while it is attached. We note that we cannot expect $\Omega_{\text{W}}$ to be the same for all pillar sizes, if the length scales of pillar, droplet and flow field, e.g. puller dipole size, are comparable (cf. Spagnolie and Lauga (2012)).

The recorded angular speeds $\omega$ are quite uniform in each orbit (see SI Fig. S7), with a slight slowdown in the second orbit, when the droplet moves on its own trail. We therefore treat $\omega$ as constant within one orbit, with $\omega\in{\omega_{1},\omega_{2}}$.

For $(\omega_{1}-\omega_{2})/\omega_{1}\ll 1$ (cf. Fig. 6 (a)), we can further assume that the orbiting period $T_{\text{orb}}$ of the droplet does not significantly increase with time, such that the droplet always experiences the same chemical gradient $\partial_{r}c(t)|_{r=R+a}$ initiated during its previous passage. We therefore approximate $\Omega_{\text{C}}$ as constant between the points where the droplet has fully crossed over onto its trail and where it leaves the pillar, and we define the detention time $t_{\text{d}}-T_{\text{orb}}$ in the second orbit as the time elapsed between these two events.

Under these approximations of constant $\Omega_{\text{W}}$, $\Omega_{\text{C}}$ and $\omega$, the detention time $t_{\text{d}}$ follows an inverse Gaussian distribution

[TABLE]

with $\lambda=\theta_{\text{e}}^{2}/2D_{\text{R}}$ and the distribution mean $\mu=\langle t_{\text{d}}\rangle=\theta_{\text{e}}/(\Omega-\omega)$.

If the deterministic torque compensates the geometry effect, the drift $\omega+\Omega$ is zero, and $\mu$ diverges. In this case, $f(t_{\text{d}})$ reverts to a Lévy distribution

[TABLE]

We have tested our approximations for the case of a straight wall without drift by calculating $f(t_{d})$ using a less approximative Fokker-Planck model Schaar et al. (2015) that includes the reflective boundary condition. The difference to the corresponding Lévy distribution is indeed too small to be resolved in experimental statistics (see SI, Sec. V).

We numerically approximate the concentration field $c(r,\alpha)$ of filled micelles in the polar coordinate system $(r,\alpha$) (Fig. 4). We assume that a droplet of radius $a$ dissolves at a constant rate and approximate its initial trail of filled micelles by a step function, $c(R<r<R+2a)=c_{0}$, which we let diffuse over time, using the micellar diffusion coefficient $D_{\mathrm{m}}\approx$100\text{,}{\mathrm{\SIUnitSymbolMicro m}}^{2}\text{,}{\mathrm{s}}^{-1}$$. We further assume the diffusion coefficients of filled and empty micelles to be similar.

Under the boundary conditions $\partial_{r}c|_{r=R}=0$ and $c(r\rightarrow+\infty)=0$ we solve the diffusion equation

[TABLE]

using a forward-marched explicit finite difference scheme. The concentration profile is computed numerically with Matlab, selecting an appropriate grid size based on the stability criteria for one-dimensional parabolic partial differential equations Grundmann (2009).

For illustration, we plot in Fig. 4(c) the rescaled filled micelle concentrations $c/c_{0}$ behind a swimmer attached to a pillar with $R=a$, diffusing from a step function at $t=0$ (dashed, blue) to a smoothed profile (solid, red) at $t=T_{\text{orb}}$. Calculated profiles and derivatives for our experimental settings are listed in the SI, Sec. IV.

For a statistical analysis of detention times, we recorded multiple sets of small numbers of monodisperse droplets of radius $a=$50\text{,}\mathrm{\SIUnitSymbolMicro m} interacting with microfluidic pillar arrays of radius $R\in\{50,75,100,250\}$\text{\,}\mathrm{\SIUnitSymbolMicro m}. Data sets for even larger pillar sizes are included in SI Sec. II.3, however, since droplets are effectively trapped at these pillars, we are not using them for detention statistics.

By restricting the analysis to “clean” pillars without residual chemorepellent from previous interactions, we retain several hundred interactions per pillar size.

We have sorted all attachment events into scattering, self-crossing and non-crossing events (Fig. 3, see SI for image analysis criteria and statistical quantities) and extracted the quantities $t_{\text{d}},v(t)$, $\omega(t)$ and $\theta(t)$. Fig. 5 contains two subsets of detention time distributions for increasing pillar sizes and the respective inverse Gaussian fits of Eq. (3). The first set (blue) contains the data range for $0<t_{\text{d}}<0.7\,T_{\text{orb}}$, during which time we expect no chemical field, $\Omega_{\text{C}}=0$. The second set (yellow, only self-crossing) contains detention times taken from the moment $T_{\text{orb}}$ where the droplet crosses its own trajectory and $\Omega_{\text{C}}=\text{const}$.

Using the extracted fitting parameters $\mu,\lambda$ from Eq. (3) and the measured $\omega$, we calculate $D_{\text{R}},\Omega_{\text{W}}$ and $\Omega_{\text{C}}$ for each pillar size, as well as the expectation value of the detention times $\langle t_{\text{d}}\rangle=\mu$, as plotted in Fig. 6.

As expected, $\langle t_{\text{d}}\rangle$ is larger in the first orbit (Fig. 6 (a)), where $\Omega_{\text{C}}=0$, and increases with increasing $R$, in fact diverging for $\Omega_{\text{C}}=0$ and $R\in\{100,250\}$\text{,}\mathrm{\SIUnitSymbolMicro m}, indicating a trapped state. The recorded mean orbital speeds $\langle|\omega|\rangle$ scale linearly with $(R+a)^{-1}$ and do not change significantly between orbits 1 and 2, justifying our approximation of an overall constant speed $v\approx$47\text{\,}\mathrm{\SIUnitSymbolMicro m}\text{\,}{\mathrm{s}}^{-1} (Fig. 6 (b)).

The extracted torques (Fig. 6 (c)) allow for a number of observations: both $\Omega_{\text{W}}$ and $-\Omega_{\text{C}}$ decrease with pillar size, the latter simply because orbiting takes longer and the chemorepulsive gradient $\partial_{r}c$ decreases by diffusion, the former possibly because the pusher type hydrodynamic interactions depend on wall curvature Spagnolie and Lauga (2012). $\Omega_{\text{C}}$ overcompensates $\Omega_{\text{W}}$ in all cases up to $R=$250\text{,}\mathrm{\SIUnitSymbolMicro m}, explaining the enhanced detachment after one orbit and the transition to fully trapped states for $R>$250\text{\,}\mathrm{\SIUnitSymbolMicro m}. Similarly, the diverging detention times in the first orbit for $R\in\{100,250\}$\text{,}\mathrm{\SIUnitSymbolMicro m}$$ are due to the fact that $\Omega_{\text{W}}$ compensates $\omega$ while $\Omega_{\text{C}}=0$.

For an estimate of the chemotactic coupling, we have plotted the fitted $\Omega_{\text{C}}$ versus the calculated filled micelle gradient, $\partial_{r}c(t)/c_{0}$, using Eq. (5) (Fig. 6 (d)). If we assume a linear dependence on the filled micelle gradient Jin et al. (2017); Taktikos et al. (2011), $\Omega_{\text{C}}=\kappa\partial_{r}c/c_{0}$, linear regression yields $\kappa\approx$130\text{,}\mathrm{rad}\text{,}\mathrm{\SIUnitSymbolMicro m}\text{,}{\mathrm{s}}^{-1} with $c_{0}$ on the order of 200 oil molecules per $\text{\,}{\mathrm{\SIUnitSymbolMicro m}}^{3}$ (see SI Sec. [II.5](#S2.SS5)), however, this should be regarded as no more than an order-of-magnitude estimate due to our limited range of accessible pillar sizes (50\text{,}\mathrm{\SIUnitSymbolMicro m}$<R<$250\text{,}\mathrm{\SIUnitSymbolMicro m}$)$ and the various approximations in our model.

In Fig. 6 (e), we plot the extracted rotational diffusion coefficients. For comparison, we have estimated the average rotational diffusion coefficient of swimmers moving between pillars by computing their mean squared rotation, $D_{\text{R}}\approx$ $0.05\text{\,}\mathrm{rad}\text{\,}{\mathrm{s}}^{-1}$ (See SI Sec. II.6). The values for droplets at pillars are quite similar, except for small pillars with chemotactic gradients, which are increased. This could be due to increased rotational noise from the backaction of the inhomogeneous chemotactic field on the surfactant coverage at the interface.

The fact that the rotational diffusion estimates for freely moving and wall attached droplets are comparable is not intuitively obvious, since variations in interfacial tension provide a major source of noise and such variations are affected by the boundary conditions of advective flow and chemical fields. Our data suggests that the differences are small enough to treat the rotational diffusion coeffient of a droplet swimmer as a meaningful physical quantity.

We have shown that from a simple curious phenomenon - droplet swimmers leaving pillars after orbiting them just once - one can independently estimate various quantities like hydrodynamic and chemorepulsive torques, diffusion coefficients and their dependence on wall curvature by a statistical analysis of detention times and their interpretation using an active Brownian swimmer model. This provides valuable insight into droplet swimmers in particular, where many of these quantities are hard to access independently or where their sensitivity to boundary conditions is still open to debate. Since the ABP model makes no assumptions about the details of the propulsion process, this assay can also be extended to probe similar quantities in various other microswimmer systems.

We acknowledge Babak V. Hokmabad for experimental, and Fabian Schwarzendahl and Babak Nasouri for theoretical advice. This project was partially funded by the MPI-DS Prandtl internship program and the MaxSynBio initiative.

I Experimental system

The microfluidic devices are fabricated with standard soft lithography techniquesQin et al. (2010). We design photomasks in a 2D AutoCad application and have them printed as a high-resolution emulsion film by an external company (128,000 dots per inch; JD Photo-Tools). From the printed photomask, we synthesize a SU-8 (Micro Resist Technology) photoresist mold on a Si wafer ((Wafer World Inc.) in a clean room environment using UV lithography. The mold is then used to generate PDMS (polydimethylsiloxane; Sylgard 184; Dow Corning) imprints. After punching in fluid inlets and outlets, the imprints are bonded to glass slides. Covalent bonding between PDMS and glass is achieved by pretreating all surfaces in an air plasma (Pico P100-8; Diener Electronic GmbH + Co. KG) for 30 s. We produce monodisperse droplets in flow-focusing microfluidic devices Thorsen et al. (2001) at high production rates ($10\text{\,}{\mathrm{s}}^{-1}30\text{\,}{\mathrm{s}}^{-1}$). Droplet radii $a$ can be tuned by flow rates and channel geometries, we used $a=$50\text{,}\mathrm{\SIUnitSymbolMicro m}$$. The syringes (Braun) are mounted on a precision microfluidic pump (NEM-B101-02B; Cetoni GmbH) and connected to the inlets and outlets with Teflon tubing (39241; Novodirect GmbH). To generate oil-in-water emulsions, the primarily hydrophobic PDMS surfaces have to be hydrophilized, first by activation via a 1:1 volumetric mixture of \ceH2O2 /HCl, then filling the channels with a silanization solution \ce(C2H2O)_nC7H18O4Si for 30 min, and finally rinsing them with milli-Q water.

For droplet production and storage we use a submicellar $0.1\,\text{~{}wt.\%}$ aqueous solution of the ionic surfactant tetradecyltrimethylammonium bromide (TTAB, critical micellar concentration $c_{\mathrm{cmc}}=0.13\,\text{~{}wt.\%}$). The concentration is sufficient to inhibit coalescence, but not high enough for oil solubilisation or droplet swimming. The oil phase is a mixture of the nematogen 4-pentyl-4’-cyano-biphenyl (5CB) and 1-Bromopentadecane (BPD) with a volume ratio of $10:1$. Pure 5CB is nematic at room temperature and deformation of the director field inside the droplet triggers helical swimming; the BPD mixture is isotropic and produces persistent Brownian swimmers.

The TTAB concentration in the experimental media is $7.5\,\text{wt}\%$, approx. $60\cdot c_{\mathrm{cmc}}$. At the beginning of each experiment, we add \approx{0.3\text{,}\mathrm{\SIUnitSymbolMicro l}} of stock oil emulsion to a $10\text{\,}\mathrm{\SIUnitSymbolMicro l}$ $7.5\,\text{wt}\%$ TTAB solution, and pipette the mixture (swimmer density ($n<$5\text{,}\mathrm{/}\mathrm{\SIUnitSymbolMicro}\mathrm{l}$$)) into the experimental container. All experiments are carried out in quasi 2D microfluidic containers. Typical PDMS cells used for experimental observations and a respective schematic drawing are shown in Fig. S1.

We observe the droplets on an Olympus IX-81 bright field microscope under 2x magnification. Images and movies are recorded by a commercial digital camera (Canon EOS 600D) at four frames per second. For one pillar size, we typically have data from more than $20$ movies recorded over $5\text{\,}\mathrm{min}10\text{\,}\mathrm{min}$ with $5$ to $20$ swimmers per cell and $1$ to $200$ individual pillar interactions per run. In the main manuscript (MS), we only included interactions at “clean” pillars, i.e. the first interaction with each pillar for each run, to avoid possible chemical contamination. Additional data for all recorded interactions are included in this supplement. Sample sizes used for the detention time statistics are listed in Table S1.

[TABLE]

II Numerical data analysis

We process video microscopy data, tracking the swimmers and extracting trajectories, using software written in-house in Python/openCV, based on a Crocker-Grier type algorithm Crocker and Grier (1996). We then identify all trajectory segments corresponding to pillar attachment (by Python script) and bin them by detention time $t_{\text{d}}$.

II.1 Criteria for trajectory segmentation

As illustrated in Fig. 2 in the MS, when a swimmer attaches or detaches from a pillar, it reorients and $\theta$ decreases or increases suddenly, leading to a peak in $|\partial_{t}\theta|$. We use this fact to numerically identify the trajectory segments corresponding to pillar attachment with the following procedure: On all recorded trajectories, we first apply distance threshold with respect to the individual pillar centres, identifying segments close to pillars. We then calculate $\alpha$ and $\theta$ along each segment (Fig. S12). The attaching and detaching points are decided by combining the thresholds $\tan(\theta)<1$ for $\theta$ and $|\tan\theta(t_{0}+0.25s)-\tan\theta(t_{0})|<0.5$ for $|\partial_{t}\theta|$, which results in an escape angle of $0.962\text{\,}\mathrm{rad}$. We have tested the robustness of this procedure by plotting and evaluating segments for multiple sample trajectories (cf. Fig. S2).

II.2 Statistical quantities

We transform the recorded Cartesian $(x,y)$ coordinates to the polar $(r,\alpha)$ system defined in Fig. 4 in the MS. Individual orbits can be identified by a path integral of $\alpha$ along the trajectory segment, detention times are given by the fixed frame rate of $4\text{\,}\mathrm{Hz}$. Detention lengths $s$ are calculated by integrating the distance between frames.

II.3 Data sets and histograms

This section contains raw histograms for a larger range of pillar sizes, $R\in\{30,50,75,100,250,\allowbreak 400,500\}$\text{,}\mathrm{\SIUnitSymbolMicro m}$$ and explains our choice of data representation in Fig. 5, MS.

We note that the data sets for $R<$50\text{,}\mathrm{\SIUnitSymbolMicro m} and $R>$250\text{\,}\mathrm{\SIUnitSymbolMicro m} represent pure scattering and trapping regimes and are therefore not usable in our analysis of chemotactically biased detachment.

We can choose between binning by detention length in the number of orbits (Fig. S3), and the detention time in seconds (Fig. S4). Binning by orbit or length clearly marks the onset of chemotactic repulsion, while binning by time is preferable if we want to fit to the time correlated stochastic noise in our Langevin model. The two approaches are only equivalent if the droplet speed is the same and absolutely constant for all swimmers; in practice, the binning by time will result in a broadening of peaks by the speed distribution, as shown in Fig. S4. To correct for this error, we have chosen to normalise all second-orbit segments to the point in time when each trajectory crosses itself. We have also excluded all interactions with pillars that have previously been hit by swimmers and might have residual chemorepellent contamination (Fig. S5).

.

II.4 Speed and angular velocity

We calculate the droplet speed $v$ via the displacement between frames, and the angular velocity $\omega$ by the change in polar angle $\alpha$, using the distance vector $\mathbf{r}$ to the pillar center. Each trajectory can be segmented into orbits of rank $n$ by $\alpha(t_{n})-\alpha(t_{n-1})=2\pi$. We have plotted distributions of $v$ and $\omega$ by pillar size and orbit rank in Fig. S6 and Fig. S7.

Since, starting from the second orbit, a droplet swims through partially consumed fuel (filled micelles), $v$ should decrease between consecutive orbits, which tendency can be seen in the respective histograms. We assume this decrease to be small in order to use a constant radial chemical gradient, however, we can use the specific angular velocities ($\omega_{1},\omega_{2}$) for our calculation of torques in the case (R=$$250\text{\,}\mathrm{\SIUnitSymbolMicro m}) where there is sufficient statistics for the second orbit.

II.5 Solubilization rate

A quantitative estimate of the droplet solubilization rate is relevant to our study in two respects: first, we have to assume a constant droplet size during each interaction, otherwise the droplet speed would change, as well as the hydrodynamic interactions. Second, we can use it to relate the chemotactic torque to the initial average concentration of oil molecules $c_{0}$ dispersed in the droplet trail.

From previous measurements under higher optical resolution (Fig. 5.1 in Krüger (2016)) we know that the change in radius over time of a 5CB droplet at room temperature in 7.5wt% TTAB is on the order of $\mathrm{d}R/\mathrm{d}t\approx$5\text{,}\mathrm{nm}\text{,}{\mathrm{s}}^{-1}$$.

Since the maximum orbital period (for R=$$250\text{\,}\mathrm{\SIUnitSymbolMicro m}) is $T_{\text{orb}}\approx$38.5\text{,}\mathrm{s}$$, the relative change in droplet size per orbit, $\frac{\Delta a}{a}(T_{\text{orb}})\approx.004$, should be negligible. Fig. S8 shows, at lower resolution and therefore quite noisy, the time evolution of droplet radii for a few trajectories from the present study. Radii decrease by less then 10% over $300\text{\,}\mathrm{s}$, which is consistent with the values given above.

To calculate the number density of dispersed 5CB oil molecules (density $1\text{\,}\mathrm{g}\text{\,}{\mathrm{cm}}^{-3}$, molar weight $250\text{\,}\mathrm{g}\text{\,}{\mathrm{mol}}^{-1}$), we calculate the volume loss per second of an $a=$50\text{,}\mathrm{\SIUnitSymbolMicro m} droplet and assume the corresponding number of molecules to disperse inside a volume given by $(2a)^{2}\cdot v\cdot$1\text{\,}\mathrm{s}.

II.6 Rotational diffusion coefficient

The rotational diffusion coefficient ($D_{\text{R}}$) can be computed in different ways. Following Howse et al. Howse et al. (2007), the mean squared displacement of an active Brownian particle is

[TABLE]

where $\tau_{R}=1/D_{\text{R}}$.

In the case of our comparatively large droplets ($a=$50\text{,}\mathrm{\SIUnitSymbolMicro m}$$), we can neglect the thermal diffusion $D_{\text{T}}$. In principle, $\tau_{r}$ can be determined by by analysing the long-time dynamics of the swimmers, calculating the MSD and identifying either the crossover time scale or fitting to diffusive limit. This is, however, not a good approach for persistent swimmers with small $D_{\text{R}}$, which within our accessible field of view ($10\text{\,}\mathrm{mm}$ $\times$ $6\text{\,}\mathrm{mm}$) do not exhibit a full diffusive crossover (Fig. S10).

We can, however, estimate $\tau_{R}$ by analyzing the decay of the angular correlation, $C(t)=\langle\mathbf{e}(\tau)\mathbf{e}\cdot(\tau+t)\rangle_{\tau}$ Saragosti et al. (2012). Inspection of the correlation decay plotted in Fig. S10 (a) for 4 trajectories shows that $\tau_{R}$ is on the order of $15\text{\,}\mathrm{s}$ to $70\text{\,}\mathrm{s}$, corresponding to a $D_{\text{R}}$ between $0.015\text{\,}\mathrm{rad}\text{\,}{\mathrm{s}}^{-1}$ and $0.067\text{\,}\mathrm{rad}\text{\,}{\mathrm{s}}^{-1}$.

We can also calculate $D_{\text{R}}$ via the mean square rotation of the swimmer by $\langle|\Delta\alpha|^{2}\rangle=2D_{\text{R}}\Delta t$ Drescher et al. (2011). In the left plot of Fig. S10, we calculate $D_{\text{r}}$ from a linear fit to $\langle|\Delta\alpha|^{2}\rangle$ for $t<$25\text{,}\mathrm{s}$$ (solid lines) for the same set of trajectories, resulting in $D_{\text{R}}$ values of $0.041$, $0.063$, $0.090$, $0.015$ respectively. Thus, we assume $D_{\text{R}}$ for droplets swimming off-wall to be on the order of $0.05\text{\,}\mathrm{rad}\text{\,}{\mathrm{s}}^{-1}$.

II.7 Data tables for the plots in Fig. 6, MS

[TABLE]

[TABLE]

III Error estimates

In Fig. 6, we provide error estimates for the quantities $\langle t_{d}\rangle$, $\omega$, $\Omega_{\text{W}}$, $\Omega_{\text{C}}$ and $D_{r}$. We directly estimate errors from experimental statistics and fitting procedures for $\theta_{\text{e}}$, $\omega$ and the fitting parameters $\lambda$ and $\mu$, and infer the values for $\langle t_{d}\rangle$, $\Omega_{\text{W}}$, $\Omega_{\text{C}}$ and $D_{r}$ via the relations $\lambda=\theta_{\text{e}}^{2}/2D_{\text{R}}$ and $\mu=\theta_{\text{e}}/(\Omega_{\text{W}}+\Omega_{\text{C}}-\omega)$ using common rules of error propagation,

[TABLE]

We note that for $R\in\{100,250\}$\text{,}\mathrm{\SIUnitSymbolMicro m}$$, $\langle t_{d}\rangle$ and the respective error diverge for the first orbit (experimentally, we cannot distinguish between inverse Gaussian and Lévy distributions). In these cases, we have not added error bars.

IV Numerical gradient estimation

To illustrate the evolution of the chemical gradient at the pillar wall, we have plotted the numerical solution of Eq. (5) for the pillar sizes discussed in the main MS in Fig. S11.

V Fokker-Planck equation for the distribution of detention times at straight walls

Here we include a Fokker-Planck solution of a swimmer interacting with a flat wall with an optional reflective boundary condition to estimate whether this condition is safe to ignore in our more approximative Langevin model.

We model our swimmer with an active Brownian particle moving at a constant speed $v$. The swimming direction is described by an overdamped Langevin equation

[TABLE]

From this equation we can a Fokker-Planck equation that accounts for the evolution in time of the one-particle probability density function, $P(\mathbf{e},t|\mathbf{e_{0}},t_{0})$, of finding an active particle, with the initial condition $P(\mathbf{e},t=t_{0}|\mathbf{e_{0}},t_{0})=\delta(\mathbf{e}-\mathbf{e_{0}})$. In the following, to lighten the notation we will use $P(\mathbf{e},t)=P(\mathbf{e},t|\mathbf{e_{0}},t_{0})$. To derive the Fokker-Planck equation, we consider the derivative of $P(\mathbf{e},t)=\langle\delta(\mathbf{e}-\mathbf{e_{0}})\rangle$ with respect to time Konotop and Vazquez (1994)

[TABLE]

The time derivative of $\langle\delta(\mathbf{e}-\mathbf{e_{0}})\rangle$ reads

[TABLE]

By using the Furutsu-Novikov-Donsker relation Sevilla and Sandoval (2015); Sevilla and Gómez Nava (2014); Kampen (1992); Klyatskin (2014); Konotop and Vazquez (1994); Frank (2005); Novikov (1965); Zinn-Justin (2002); Vachier and Mazza (2019)

[TABLE]

the Fokker-Planck equation reads

[TABLE]

where $\text{L}_{e}^{2}$ is the angular part of the Laplacian. Let us take a step back and assume the drift $\Omega$ be null. Then we can rewrite Eq. (S1) as

[TABLE]

This system can be then mapped onto the quantum mechanical problem of the rigid rotator and by using a path integral representation Zinn-Justin (2002) it is possible to find the solution of the Fokker-Planck equation $P(\mathbf{e},t)$ in $d$ dimensions. By using the bracket notation $P(\mathbf{e},t)$ reads

[TABLE]

From Zinn-Justin (2002), the path integral representation for the probability distribution $P(\mathbf{e},t)$ in a free 2D plane can be found and reads

[TABLE]

In order to calculate the distribution of detention times, we need to impose an absorbing boundary condition at the escape angle $\theta_{\text{e}}$ , $[P(\theta_{\text{e}},t)=0]$. By using the method of images, the probability distribution can be found and reads

[TABLE]

The survival probability is given by

[TABLE]

The detention time follows a Lévy distribution

[TABLE]

We now apply a reflective boundary condition $[\partial_{\theta}P(\theta,t)|_{\theta_{\text{r}}}=0]$. By using the method of images, the survival probability is given by

[TABLE]

The detention time is given by $f(t)=-\frac{\partial}{\partial t}S(t)$ and reads

[TABLE]

We plot Eq. (S2) and Eq. (S3) with parameters relevant to our experimental system ($\theta_{0}=0$, $\theta_{\text{e}}=-0.962$, $\theta_{\text{r}}=\pi$, and $D_{\text{R}}=0.05\,\text{rad}^{2}/\text{s}$) in Fig. S13; the difference between two curves is indeed negligible.

We note that for different $D_{\text{R}}$, the distribution curves will be similar in shape, as Eq. (S2) can be written as

[TABLE]

with $\bar{t}=D_{\text{R}}t$. This also holds for Eq. (S3).

VI Detention time distribution at large pillars

The detention time distributions at large pillars have a large mean, such that $f(t_{d})$ is close to the Lévy distribution for a flat wall, as derived above. In Fig. S14 we show the normalised detention time histograms at large pillars (radii $400\text{\,}\mathrm{\SIUnitSymbolMicro m}$ and $500\text{\,}\mathrm{\SIUnitSymbolMicro m}$), and fits with both inverse Gaussian and Lévy distributions. The differences are within the experimental error.

VII Supporting movie