Swimming of microorganisms in quasi-2D membranes

Microscopic biological organisms have adapted to a viscous world in which their inertia is inconsequential to their locomotion. For a typical micro-scale organism the Reynolds number $\Rey=\rho UL/\eta$ is small. For example, Escherichia coli has a characteristic length $L\sim 10$ $\mu$m and a characteristic swimming speed $U\sim 10\,\mu$m/s in water (density $\rho\approx 10^{3}$ $\rm{kg/m^{3}}$ and dynamic viscosity $\eta\approx 10^{-3}$ Pa s), leading to a negligibly small Reynolds number $\Rey=\rho UL/\eta\sim 10^{-5}$–$10^{-4}$ (Purcell, 1977). At this scale, a swimmer reacts instantaneously to any forces, oblivious to any history of prior dynamics (Purcell, 1977; Childress, 1981; Lauga & Powers, 2009), and so the primary method of locomotion for macroscopic swimmers such as fish and humans, which relies on trading momentum with the fluid to generate a propulsive force, is ineffective here. Rather, the net translations of a microorganism are determined by the sequence of configurations it adopts to swim, independent of its deformation rate. Microswimmers must continually paddle or deform their bodies in a swimming pattern with non-reciprocal forward and reverse strokes to manipulate the drag forces for propulsion (Purcell, 1977; Childress, 1981; Lauga & Powers, 2009).

It is common for microorganisms to swim in geometrically confined media: in channels, near surfaces and interfaces, and in films. A significant amount of theoretical and experimental work has been devoted to studying the effects of confinement on the motion of microscopic swimmers near solid walls (Pedley & Kessler, 1987; Lauga et al., 2006; Berke et al., 2008; Drescher et al., 2009; Li & Tang, 2009; Or & Murray, 2009; Or et al., 2011; Crowdy & Or, 2010; Li et al., 2011; Spagnolie & Lauga, 2012; Molaei et al., 2014; Ishimoto et al., 2016), near fluid-fluid interfaces (Guasto et al., 2010; Di Leonardo et al., 2011; Wang & Ardekani, 2013; Lopez & Lauga, 2014; Masoud & Stone, 2014; Stone & Masoud, 2015), and in thin fluid layers atop a solid substrate (Lambert et al., 2013; Mathijssen et al., 2016a, b; Ota et al., 2018).

Motivated by recent experimental and theoretical studies of bacteria swimming in biofilms, in freely-suspended thin films (Aranson et al., 2007; Sokolov et al., 2007), and on active proteins mimicking biological swimmers in lipid membranes (Huang et al., 2012), we study here the hydrodynamics of swimming microorganisms in a thin membrane. We treat the membrane as a continuous incompressible viscous fluid film of very small thickness. Flow fields in such a membrane are uniform throughout the thickness of the membrane. In contrast to a thin-film model that involves integration of three-dimensional (3D) hydrodynamic equations across the thickness of the film, our membrane model is intrinsically two-dimensional (2D), in the sense that the motion of molecules within the membrane in the direction normal to the plane of the membrane is forbidden. The membrane is embedded in a 3D fluid of different viscosity. The motion of a swimmer in the membrane generates flows both in the membrane and in the surrounding fluid (see figure 1). While there have been some recent investigations of the hydrodynamics of swimmers in a thin layer of fluid sandwiched between fluids of a different viscosity (Leoni & Liverpool, 2010; Rower et al., 2019), this problem is still largely unexplored.

As was demonstrated by Saffman & Delbrück (1975) and Saffman (1976), the amount of momentum imparted by the membrane to the bulk fluid is controlled by a hydrodynamic length scale, the so-called Saffman length $\ell_{S}$, given by the ratio of the 2D membrane viscosity $\eta_{m}$ to the 3D viscosity of the bulk fluid $\eta$, $\ell_{S}=\eta_{m}/(2\eta)$. If the membrane is perturbed by a localized force (applied in the plane of the membrane) at a point $\boldsymbol{x}$, for distances (measured from $\boldsymbol{x}$) much smaller than the Saffman length, the effect of the flows in the bulk on the membrane hydrodynamics is negligible. In this region the fluid velocity in the membrane decays slowly (logarithmically) with distance, as in purely 2D fluids. On the other hand, for distances much larger than the Saffman length, the contribution of the bulk fluid to the membrane dynamics is significant, and the membrane flow field decays inversely with the distance, a behavior consistent with 3D dynamics.

Levine & MacKintosh (2002) (LM) derived a Green function for a more general case of viscoelastic membranes. In the case of a purely viscous membrane that we consider here, there is no elastic response of the membrane, and a disturbance caused by a force results in the velocity field alone. The velocity of the membrane at position $\boldsymbol{x}^{\prime}$ due to an in-plane, localized force $\boldsymbol{f}(\boldsymbol{x})=\boldsymbol{f}\delta(\boldsymbol{x})$ is determined by the LM response tensor $\mbox{\boldmath$\alpha$}({\boldsymbol{x}}-{\boldsymbol{x}}^{\prime})$,

$${\boldsymbol{v}}({\boldsymbol{x}})=\frac{1}{4\pi\eta_{m}}\mbox{\boldmath$\alpha$}({\boldsymbol{x}}-{\boldsymbol{x}}^{\prime})\bcdot{\boldsymbol{f}}({\boldsymbol{x}}^{\prime}).$$

(1)

Here ${\boldsymbol{x}}$ and ${\boldsymbol{x}}^{\prime}$ are in-plane vectors with components $(x,y)$ and $(x^{\prime},y^{\prime})$, respectively (refer also to figure 1 for our choice of the coordinate system). The response function $\mbox{\boldmath$\alpha$}({\boldsymbol{x}}-{\boldsymbol{x}}^{\prime})$ in equation (1) plays the role of the Oseen tensor in $3$D hydrodynamics.

As was shown by LM, the response function may be split into ‘parallel’ and ‘transverse’ contributions. In the component form we have

$$\alpha_{\alpha\beta}({\boldsymbol{x}})=\alpha_{\|}(|{\boldsymbol{x}}|)\hat{x}_{\alpha}\hat{x}_{\beta}+\alpha_{\perp}(|{\boldsymbol{x}}|)[\delta_{\alpha\beta}-\hat{x}_{\alpha}\hat{x}_{\beta}],$$

(2)

where $\alpha,\beta=x,y$, and $\hat{x}_{\alpha}$ is the $\alpha$ component of the unit vector $\hat{{\boldsymbol{x}}}={\boldsymbol{x}}/|{\boldsymbol{x}}|$. In our notation, $\alpha_{\alpha\beta}$ corresponds to $-i\omega\alpha_{\alpha\beta}$ in the LM theory. The scalar functions $\alpha_{\|}$ and $\alpha_{\perp}$ are given by

	$\displaystyle\alpha_{\|}(\kappa)$	$\displaystyle=$	$\displaystyle\frac{\pi}{\kappa}{\boldsymbol{H}}_{1}(\kappa)-\frac{2}{\kappa^{2}}-\frac{\pi}{2}[Y_{0}(\kappa)+Y_{2}(\kappa)]$		(3)
	$\displaystyle\alpha_{\perp}(\kappa)$	$\displaystyle=$	$\displaystyle\pi{\boldsymbol{H}}_{0}(\kappa)-\frac{\pi}{\kappa}{\boldsymbol{H}}_{1}(\kappa)+\frac{2}{\kappa^{2}}-\frac{\pi}{2}[Y_{0}(\kappa)-Y_{2}(\kappa)],$

where ${\boldsymbol{H}}_{\nu}$ are Struve functions and $Y_{\nu}$ are Bessel functions of the second kind (Abramowitz & Stegun, 1965); $\kappa=|{\boldsymbol{x}}|/\ell_{S}$ is the non-dimensionalized distance between the point of application of the force and the point where the membrane velocity response is measured. Both $\alpha_{\|}(\kappa)$ and $\alpha_{\perp}(\kappa)$ diverge logarithmically as $\kappa\to 0$, while for large $\kappa$ we have $\alpha_{\|}(\kappa)\sim 1/\kappa$ and $\alpha_{\perp}(\kappa)\sim 1/\kappa^{2}.$

In the small Reynolds number regime the inertia term in the Navier-Stokes equation can be neglected. We assume that the membrane has thickness $h$ and choose a coordinate system with the origin at the top surface of the membrane with $z=0$ (therefore, the bottom side of the membrane is at $z=-h$). The dynamics of the membrane embedded in a bulk fluid is governed by a modified Stokes equation and the incompressibility condition (Saffman, 1976),

$\displaystyle-\bnabla p+\frac{\eta_{m}}{h}\bnabla^{2}\boldsymbol{v}+\frac{2\boldsymbol{f}}{h}=0,\qquad\bnabla\cdot\boldsymbol{v}=0,$

(4)

where $p$ and $\boldsymbol{v}$ are the pressure and velocity fields of the membrane. The flows in the membrane set the bulk fluid into motion. The resulting flows in the embedding fluid, in their turn, exert traction on the membrane. Since the membrane is only a few molecular layers thick, the traction due to the bulk fluid produces a flow that is uniform throughout the thickness of the membrane, i.e. the fluid velocity does not depend on $z$ for $-h<z<0$. In equation (4) the coupling between the membrane and the bulk fluid is described by the force per unit volume $2\boldsymbol{f}/h$, with

$$\boldsymbol{f}=\eta\frac{\partial\boldsymbol{v}^{(3D)}}{\partial z}{\Big{|}}_{z=0},$$

(5)

where $\boldsymbol{v}^{(3D)}$ is the bulk fluid velocity. The factor of two in equation (4) is due to an equal force $\boldsymbol{f}$ acting on the bottom side of the membrane. Equation (4) can be written in a more compact form that we will use later,

$$\bnabla\cdot\boldsymbol{\sigma}=-\frac{2\boldsymbol{f}}{h},$$

(6)

where $\boldsymbol{\sigma}$ is the stress field of the membrane.

Being inertialess and swimming in the absence of external influences, a swimmer must maintain a zero net force $\boldsymbol{F}(t)$ and a zero torque $\boldsymbol{L}(t)$ on its body at every time instant,

	$\displaystyle\boldsymbol{F}(t)$	$\displaystyle=$	$\displaystyle\int_{S}\boldsymbol{\sigma}\cdot\boldsymbol{n}\,\mathrm{d}S,$		(7)
	$\displaystyle\boldsymbol{L}(t)$	$\displaystyle=$	$\displaystyle\int_{S}\boldsymbol{x}\times(\boldsymbol{\sigma}\cdot\boldsymbol{n})\,\mathrm{d}S,$		(8)

where the integration is over the surface of the swimmer and $\boldsymbol{n}$ is a unit vector normal to the surface and pointing away from the swimmer.

Many microorganisms such as spermatozoa, E. Coli, and Caulobacter crescentus swim by moving thin extensions (flagella) on their bodies. Some critters, like Paramecium, are covered in thousands of short hair-like appendages called cilia and propel themselves through a coordinated beating of these cilia. Since our primary goal is to study how the confinement to the plane of a membrane affects the swimming dynamics of a microorganism (rather than a detailed study of a particular microorganism), we consider here only minimal theoretical models of flagellated and ciliated microorganisms.

In §2 we consider a headless, infinitely long ‘flagellum’ of infinitesimally small thickness propagating planar sinusoidal waves along its body. This is a one-dimensional analog of the Taylor swimmer (Taylor, 1951), an infinite plane in viscous fluid passing transverse sinusoidal waves. We recover Taylor’s result for the swimming velocity in the limiting case of a pure 2D hydrodynamics (the membrane in vacuum). We find that the membrane incompressibility condition imposes a constraint on the fluid dynamics that allows the flagellum to achieve much higher swimming speed than in pure 2D and 3D fluids for large ratios of the wavelength to the Saffman length.

In §3 we study the propulsion of a flagellum of finite length and find its swimming speed and efficiency. In §2 and §3 we apply the boundary-element method (BEM) that two of us have recently developed in work on hydrodynamic interaction of inclusions in freely-suspended smectic films (Qi et al., 2014; Kuriabova et al., 2016; Qi et al., 2017).

In §4 we formulate the Lorentz reciprocal theorem for a quasi-2D fluid and derive an equation for the swimming speed. We discuss the advantage of the method based on the Lorentz reciprocal theorem over the BEM for studying microscopic swimmers with swimming patterns that do not change the overall shape of the swimmers’ bodies (for example, swimmers propagating longitudinal compressive waves along their bodies). As an example of such a critter, we consider a simple model of a two-dimensional ‘squirmer’, a disk with a prescribed tangential velocity along its circumference. We find that unlike its flagellated counterpart, a squirmer does not benefit from the presence of the bulk fluid: its swimming speed is lower than that in a purely 2D fluid.

In §5 we discuss our results and suggest further directions of investigation.

The Taylor swimmer (Taylor, 1951) is an infinite swimming 2D sheet in a 3D viscous fluid that propagates transverse waves of amplitude $b$ and wave speed $c=\omega/q$. In a frame moving with the swimming sheet (co-moving frame) the shape of the sheet is described by

$$y=b\sin(qx-\omega t)=b\sin\xi,$$

(9)

where $\xi=qx-\omega t$ denotes the wave phase.

Taylor showed that the sheet with such a one-dimensional modulation travels, relative to the fluid at infinity, with a speed $U/c=(bq)^{2}/2+O((bq)^{4})$ in the direction opposite to the wave velocity. We consider here a one-dimensional analog of the Taylor swimmer: an infinitely long, infinitesimally thin flagellum confined to the plane of a viscous membrane embedded in bulk fluid (see figure 1) with prescribed motion given by equation (9), with $x$ and $y$ parametrizing the shape of the flagellum.

A swimming velocity of an infinitely long flagellum is time-independent. Indeed, two snapshots of the waving flagellum taken at the same point on the $x$-axis differ only by a shift $\Delta x$ along the $x$-axis. Thus, a temporal shift at a fixed point $x$ is equivalent to a spatial displacement along the flagellum. A swimmer moving as a whole has the same translational velocity along its entire length. The swimming velocity, being invariant with respect to translations along the $x$-axis must be invariant with respect to translations in time as well. Therefore, we can calculate the swimming speed for a single time instant and set $t=0$ in equation (9).

As in (Taylor, 1951) we consider the case of an inextensible flagellum. In order to calculate the portion of the swimmer’s velocity due to its distortion, we calculate the position of a material point of the flagellum as a function of time. In a frame moving with the wave (with speed $c$ relative to the co-moving frame) the shape of the flagellum does not change. The Cartesian coordinates for this frame are $(x^{\prime},y)$, where $x^{\prime}=x-ct$. In this reference frame a material particle of the flagellum travels a distance $\Lambda$ equal to the arclength of the flagellum spanned by one (linear) wavelength $\lambda$ during one period of oscillation $T=2\pi/\omega$,

$$\Lambda=\frac{1}{q}\int_{0}^{2\pi}\sqrt{1+(bq)^{2}\cos^{2}\xi}\mathrm{d}\xi.$$

(10)

We will call $\Lambda$ the arcwise wavelength. The material particle’s speed is, therefore,

	$\displaystyle C$	$\displaystyle=$	$\displaystyle\frac{c}{2\pi}\int_{0}^{2\pi}\sqrt{1+(bq)^{2}\cos^{2}\xi}\mathrm{d}\xi$		(11)
		$\displaystyle\approx$	$\displaystyle c\left(1+\frac{1}{4}b^{2}q^{2}-\frac{3}{64}b^{4}q^{4}\right).$		(12)

To determine the position of a material particle of the flagelleum, we define the material coordinate $S$ to be the arclength coordinate $s$ of a material point at $t=0$. In the frame in which the nodes of the wave are fixed, arclength is related to the Cartesian coordinate $x^{\prime}$ by

	$\displaystyle s$	$\displaystyle=$	$\displaystyle\int_{0}^{x^{\prime}}\sqrt{1+(bq)^{2}\cos^{2}(qx^{\prime})}\mathrm{d}x^{\prime}$		(13)
		$\displaystyle\approx$	$\displaystyle x^{\prime}+\frac{b^{2}q}{8}\left[2qx^{\prime}+\sin(2qx^{\prime})\right]-\frac{b^{4}q^{3}}{256}\left[12qx^{\prime}+8\sin(2qx^{\prime})+\sin(4qx^{\prime})\right].$		(14)

Reverting the series leads to

$$x^{\prime}\approx s-\frac{b^{2}q}{8}\left[2qs+\sin(2qs)\right]+\frac{b^{4}q^{3}}{256}\left[28qs+16qs\cos(2qs)+16\sin(2qs)+5\sin(4qs)\right].$$

(15)

Using $y=b\sin qx^{\prime}$ and $s=S-Ct$ leads to the position of the material point labeled by $S$ as a function of time $t$:

	$\displaystyle x$	$\displaystyle\approx$	$\displaystyle S-\frac{b^{2}q}{8}\left[2qS+\sin 2(qS-\omega t)\right]$		(16)
			$\displaystyle+\frac{b^{4}q^{3}}{256}\left[28qS+16qS\cos 2(qS-\omega t)+16\sin 2(qS-\omega t)+5\sin 4(qS-\omega t)\right],$
	$\displaystyle y$	$\displaystyle\approx$	$\displaystyle b\sin(qS-\omega t)$		(17)
			$\displaystyle-\frac{b^{3}}{16}\left[4q^{3}S\cos(qS-\omega t)+q^{2}\sin(qS-\omega t)+q^{2}\sin 3(qS-\omega t)\right].$

In the co-moving frame the components of a material particle’s velocity $\boldsymbol{u}_{S}$ are given by

$$\boldsymbol{u}_{S}=\left(\left.\frac{\partial{x}}{\partial t}\right|_{S},\left.\frac{\partial{y}}{\partial t}\right|_{S}\right).$$

(18)

For an arbitrary value of $b$, the material particle’s velocity can be calculated numerically using

	$\displaystyle u_{S,x}$	$\displaystyle=$	$\displaystyle-C\cos\theta_{S}+c$		(19)
	$\displaystyle u_{S,y}$	$\displaystyle=$	$\displaystyle-C\sin\theta_{S},$		(20)

with $\displaystyle\tan\theta_{S}=\left.\frac{\partial y}{\partial x}\right|_{S}$. The total velocity of a material particle relative to the fluid at infinity (for $y\to\pm\infty$) is the sum of the surface disturbance and swimming velocities, $\boldsymbol{u}_{S}+\boldsymbol{U}$.

The linearity of Stokes equations allows us to model the fluid velocity field in the membrane as a superposition of fluid velocities due to a (yet unknown) force density $\boldsymbol{f}(\boldsymbol{x})$ along the flagellum:

$$\boldsymbol{v}(\boldsymbol{x})=\frac{1}{4\pi\eta_{m}}\int_{\Gamma}\mbox{\boldmath$\alpha$}(\boldsymbol{x}-\boldsymbol{x}^{\prime})\bcdot\boldsymbol{f}(\boldsymbol{x}^{\prime})\mathrm{d}\boldsymbol{x}^{\prime},$$

(21)

where the integration is along the (infinite) contour of the flagellum.

The spatial periodicity of the flagellum modulation implies the invariance of the flow field and the force density $\boldsymbol{f}(\boldsymbol{x})$ in equation (21) under translations along the $x$-axis by an integer multiple of the wavelength $\lambda=2\pi/q$. Defining $x_{m}=x+m\lambda$, for all integers $m$, the integration on the RHS of equation (21) can be reduced to integration over a single wavelength. We impose a no-slip boundary condition on the surface of the flagellum by setting the fluid velocity equal to the velocity of the material point on the flagellum,

$$\boldsymbol{u}_{S}(\boldsymbol{x})+\boldsymbol{U}=\frac{1}{4\pi\eta_{m}}\int_{\Gamma_{0}}\sum_{m=-\infty}^{\infty}\mbox{\boldmath$\alpha$}(\boldsymbol{x}-\boldsymbol{x}_{m}^{\prime})\bcdot\boldsymbol{f}(\boldsymbol{x}^{\prime})\mathrm{d}\boldsymbol{x}^{\prime},$$

(22)

where $\boldsymbol{x}$ and $\boldsymbol{x}^{\prime}$ indicate the points on the flagellum that belong to a one-wavelength ‘window’ $\Gamma_{0}$ and $\boldsymbol{x}_{m}^{\prime}=(x^{\prime}+m\lambda,y^{\prime})$. The surface disturbance velocity $\boldsymbol{u}_{S}(\boldsymbol{x})$ on the LHS of equation (22) is given by equation (18).

To close the system of equations for the force density $\boldsymbol{f}(\boldsymbol{x})$ and the swimming velocity $\boldsymbol{U}$, we also require the net force on the flagellum be equal to zero,

$$\int_{\Gamma_{0}}\boldsymbol{f}({\boldsymbol{x}})\mathrm{d}\boldsymbol{x}=0.$$

(23)

We solved equations (22) and (23) numerically in Matlab by splitting the integration path into $N$ straight-line segments of equal length $\Delta s$ and replacing the line integrals in equations (22) and (23) by summation over the segments,

$$\boldsymbol{u}_{S}(\boldsymbol{x}_{i})+\boldsymbol{U}=\frac{1}{4\pi\eta_{m}}\sum_{j=1}^{N}\sum_{m=-M}^{M}\mbox{\boldmath$\alpha$}(\boldsymbol{x}_{i}-\boldsymbol{x}_{mj})\bcdot\boldsymbol{f}(\boldsymbol{x}_{j})\Delta s,$$

(24)

$$\sum_{j=1}^{N}\boldsymbol{f}(\boldsymbol{x}_{j})=0.$$

(25)

In equations (24) and (25) $\boldsymbol{x}_{i(j)}$ are the coordinates of the segments’ midpoints. In equation (24) we introduced the truncation parameter $M$ for the infinite sum over the wavelengths.

For a starting value of parameter $M$ (usually $M=10$), we ran computations for five different values of parameter $N$ in the range $300$–$1000$ and then extrapolated our results for the swimming velocity to $\Delta s\to 0$ ($N\to\infty$). We then gradually increased the value of $M$ and repeated the computations until the solution for $\boldsymbol{U}$ converged, showing changes smaller than $0.5\%$ with further increase of the number of terms in the sum over $m$. The computations required increasingly more terms in the sum over $m$ for large amplitudes ($bq>1$) and large Saffman lengths ($\lambda/\ell_{S}\ll 1$) due to strong long-range hydrodynamic interactions between segments of the flagellum in this (nearly 2D) regime, and correspondingly large contribution to the flow field by the forces $\boldsymbol{f}(\boldsymbol{x}_{j})$ separated by multiple wavelengths along the flagellum.

We paid special attention to the diagonal term with $m=0$ and $i=j$ in equation (24). This term gives the fluid velocity in the close proximity of a localized force $\boldsymbol{f}(\boldsymbol{x})$. The response tensor $\mbox{\boldmath$\alpha$}(\boldsymbol{x})$ diverges due to logarithmic singularities in the functions $\alpha_{\|}(\boldsymbol{x})$ and $\alpha_{\perp}(\boldsymbol{x})$ in the limit $\boldsymbol{x}\to 0$ [see equations (3) and (1)]. In the close proximity of a localized force the fluid velocity is parallel to the force, and is, therefore, determined by the parallel component of the response function $\alpha_{\|}(\boldsymbol{x})$. We expanded $\alpha_{\|}(\boldsymbol{x})$ about $\boldsymbol{x}=0$ and performed integration analytically over $\Delta s$ in the vicinity of $\boldsymbol{x}=0$. Therefore, for the diagonal term on the RHS of equation (24), which we denote as $\boldsymbol{{\cal A}}^{m=0}_{i=j}$, we have

	$\displaystyle\boldsymbol{{\cal A}}^{m=0}_{i=j}(\boldsymbol{x}_{i})$	$\displaystyle=$	$\displaystyle\frac{1}{4\pi\eta_{m}}\boldsymbol{f}(\boldsymbol{x}_{i})\int_{-\Delta s/2}^{\Delta s/2}\alpha_{\|}(z)\mathrm{d}z$		(26)
		$\displaystyle=$	$\displaystyle\frac{1}{4\pi\eta_{m}}\boldsymbol{f}(\boldsymbol{x}_{i})\,2\lim_{\varepsilon\to 0}\int_{\varepsilon}^{\Delta s/2}\left[\frac{1}{2}-\gamma+\frac{2z}{3\ell_{S}}+\log\frac{2\ell_{S}}{z}\right]\mathrm{d}z$
		$\displaystyle=$	$\displaystyle\frac{1}{4\pi\eta_{m}}\boldsymbol{f}(\boldsymbol{x}_{i})\,\Delta s\left[\frac{3}{2}+\frac{\Delta s}{6\ell_{S}}-\gamma+\log\left(\frac{4\ell_{S}}{\Delta s}\right)\right],$

where $\gamma=0.5772$ is the Euler constant.

Our computations confirm that an infinitely long flagellum has a non-vanishing component of the swimming velocity only along the $x$-axis, as expected by symmetry. In figure 2 we plot the swimming speed as a function of the dimensionless amplitude $bq$ for a range of wavelengths scaled by the Saffman length. In the limit of pure 2D hydrodynamics, which corresponds to large Saffman lengths (and small scaled wavelengths, $\lambda/\ell_{S}\ll 1$), the energy dissipation occurs primarily in the membrane, and the membrane’s viscous drag on the flagellum makes the main contribution to the flagellum’s propulsion. In this limit, our computations are in good agreement with the 2D problem of the Taylor swimming sheet, as expected because a Taylor ‘string’ in a thin very viscous membrane is equivalent to the Taylor sheet in 3D bulk fluid. For small amplitudes ($bq\ll 1$), we recover Taylor’s leading order perturbative solution $U/c=(1/2)(bq)^{2}$. For larger amplitudes $bq$ (and $\lambda/\ell_{S}\ll 1$), our calculated swimming speed is in agreement with recent analytic and computational results of Sauzade and coworkers (Sauzade et al., 2011).

In figure 2 the black dashed curve corresponds to a flagellum swimming in a pure 2D membrane (no bulk fluid surrounding the membrane). We briefly outline our computations for this limiting case in Appendix A. The computations are similar to the boundary integral approach demonstrated in (Sauzade et al., 2011) with the only difference that we neglected the double layer potential contribution. The solid black curve in figure 2 corresponds to the local drag theory for an infinitely long flagellum passing waves of moderately large amplitudes in a 3D unbounded fluid (Gray & Hancock, 1955).

As can be seen in figure 2, our BEM computations predict that for wavelengths larger than the Saffman length ($\lambda/\ell_{S}>1$) the swimming speed in a quasi-2D membrane exceeds that in purely 2D and 3D fluids. For qualitative explanation of this result we compare our BEM computations with the local drag model of Gray & Hancock (1955). In the local drag approximation, one assumes that the viscous drag force on a small segment of the flagellum is proportional to the segment’s velocity, and the total drag on the swimmer is a sum over these local drag forces. Thus, the local drag approximation does not explicitly take into account the long-range hydrodynamic interactions between distant segments of the flagellum. The local drag forces for the motion of a rod-like segment parallel and perpendicular to its geometrical axis are given by $F_{\|}=\zeta_{\|}v_{\|}$ and $F_{\perp}=\zeta_{\perp}v_{\perp}$, respectively, with the drag coefficients $\zeta_{\|}$ and $\zeta_{\perp}$.

We expect our BEM results to be in qualitative agreement with the local drag approximation in the limit of $bq\ll 1$ and $\lambda/\ell_{S}\gg 1$. For small amplitudes $bq$ the segments of an inextensible flagellum separated by large contour distances ($>\lambda$) do not come too close to each other, and for the wavelengths larger than the Saffman length the spatial decay of the flow field is faster ($\sim 1/r$), in comparison with slower (logarithmic) decay rate for $\lambda/\ell_{S}\ll 1$. Thus in the regime of $bq\ll 1$ and $\lambda/\ell_{S}\gg 1$, the cooperativity effect between distant segments of the flagellum is expected to be small. Gray & Hancock (1955) obtained the swimming velocity of an infinitely long and thin flagellum to the leading order of amplitude $bq$,

$$\frac{U}{c}=\frac{(bq)^{2}}{2}\left(\frac{\zeta_{\perp}}{\zeta_{\|}}-1\right).$$

(27)

In 3D, the ratio of the drag coefficients for an infinitely thin rod is $\zeta_{\perp}/\zeta_{\|}=2$. For inclusions in quasi-2D membranes, the ratio $\zeta_{\perp}/\zeta_{\|}$ depends on the Saffman length. In the inset of figure 2 we plot our BEM results for $(\zeta_{\perp}/\zeta_{\|})_{\rm estim}\equiv 1+2U/(cb^{2}q^{2})$ as a function of $\lambda/\ell_{S}$. According to equation 27, $(\zeta_{\perp}/\zeta_{\|})_{\rm estim}$ should give us an estimate for the local drag anisotropy. For the plot in the inset we chose a small amplitude $bq=10^{-3}$, when the comparison with the local drag calculation of Gray and Hancock is justified. As can be seen in the inset of figure 2, the effective ratio $(\zeta_{\perp}/\zeta_{\|})_{\rm estim}$ grows logarithmically with $\lambda/\ell_{S}$ for $\lambda/\ell_{S}\gg 1$.

This result is in qualitative agreement with the work of Levine and collaborators (Levine et al., 2004), some of which we summarize here. Levine et al. studied the drag coefficients for a rod-like inclusion of length $L$ moving in a quasi-2D membrane, and showed that for rod-like inclusions of lengths smaller than the Saffman length ($L/\ell_{S}\ll 1$), where the viscous dissipation occurs primarily in the membrane, the dependence of the drag coefficients on the size and orientation of the rod is weak: $\zeta_{\perp}/\zeta_{\|}\to 1$. For longer rods with $L/\ell_{S}\gg 1$, the dissipation is governed by the 3D fluid surrounding the membrane, and the drag coefficients show a stronger dependence on the size of the rods. Levine et al. found that the drag coefficient $\zeta_{\|}$ for a thin rod in a quasi-2D membrane is qualitatively similar to that in 3D and is given by

$$\zeta_{\|}=\frac{2\pi\eta L}{\ln\left(0.43L/\ell_{S}\right)}.$$

(28)

However, the dependence of $\zeta_{\perp}$ on $L$ in a quasi-2D membrane is very different from that in 3D:

$$\zeta_{\perp}=2\pi\eta L;$$

(29)

$\zeta_{\perp}$ depends on $L$ linearly, without the logarithmic factor in the denominator. The linear dependence of $\zeta_{\perp}$ on $L$ indicates the local character of the drag and the effective absence of hydrodynamic interactions between different sections of the rod.

As emphasized by Levine et al. (2004), this behavior of $\zeta_{\perp}$ arises from the incompressibility of the membrane, $\bnabla_{\perp}\bcdot\boldsymbol{v}_{\perp}=0$, where the symbol $\perp$ denotes differentiation ‘in-plane’ and the components of the velocity field in the plane. In the case of a rod moving perpendicular to its long axis in a 3D fluid, the fluid can flow past the rod by moving over and under it. In this flow pattern, the in-plane part of the incompressibility condition does not vanish: $\partial_{x}v_{x}+\partial_{y}v_{y}\neq 0$. A 2D version of such a flow (its projection on the $xy$-plane) is impossible due to the membrane incompressibility. In a quasi-2D membrane, the fluid moves the long way around the rod, and the flow extends over distances comparable to the largest rod dimension $L$. The membrane incompressibility constraint only affects the perpendicular drag on a filament. A segment of filament being dragged parallel to its long axis does not produce divergent flows in a simple fluid, and thus the flow character is unchanged by the presence of a membrane.

Therefore, for long wavelengths ($\lambda/\ell_{S}\gg 1$), the membrane incompressibility is expected to lead to a logarithmic growth of the flagellum’s effective drag anisotropy, $\zeta_{\perp}/\zeta_{\|}\propto\log(\lambda/\ell_{S})$. An organism that relies on the drag anisotropy for propulsion would achieve greater swimming speeds in a quasi-2D membrane than in pure 2D or 3D fluids, as is confirmed by our BEM computations.

We also applied the BEM to the case of an inextensible, headless, infinitely thin flagellum of finite length. As in the case of an infinitely long flagellum, the motion of the swimmer is prescribed by a sinusoidal modulation, $y(s,t)=b\sin(qx(s)-\omega t+\phi_{0})$. Here $s$ is the arc length along the flagellum measured from the flagellum’s hypothetical ‘head’ and $\phi_{0}$ is the initial phase. At every time instant, the shape of the flagellum is described by the curve $\boldsymbol{X}(s,t)$, where $\boldsymbol{X}(s,t)=(X(s),Y(s))=(x(s),y(s,t)-y(0,t))$ (see figure 3). The unit tangent to the curve is $\boldsymbol{T}(s)=(\mathrm{d}X/\mathrm{d}s,\mathrm{d}Y/\mathrm{d}s)$.

In the frame of the flagellum, a material point at position $s$ moves with velocity $\boldsymbol{u}_{S}(s,t)=\partial\boldsymbol{X}(s,t)/\partial t$. As was demonstrated by Higdon (1979), in the case of transverse waves propagating along the flagellum, the material particle’s velocity can be calculated in a different manner. In a reference frame moving with the wave, the shape of the flagellum is given by $\boldsymbol{X}^{w}(s-Ct)$, where $C$ is the arcwise speed that we introduced in equation (11). Following Higdon, we note that

$$X^{w}(s+\Lambda)=X^{w}(s)+\lambda,\qquad Y^{w}(s+\Lambda)=Y^{w}(s),$$

(30)

where $\Lambda$ is the arcwise wavelength (see equation (10)) and $\lambda$ is the linear wavelength. The tangential vectors are identical in the body and wave frames since the wave frame simply translates with respect to the body frame and does not undergo rotation. Thus, $\boldsymbol{T}^{w}(s-Ct)=\boldsymbol{T}(s-Ct)$. The velocity of a material particle at $s$ in the wave frame is calculated as $\boldsymbol{u}^{w}(s,t)=\partial\boldsymbol{X}^{w}(s-Ct)/\partial t=-C\partial\boldsymbol{X}^{w}/\partial s=-Q\boldsymbol{T}(s-Ct)$. The velocity of the ‘head’ in the wave frame is $-C\boldsymbol{T}(s-Ct)\Big{|}_{s=0}=-C\boldsymbol{T}(-Ct)$. Therefore, the wave frame translates with respect to the body frame with velocity $C\boldsymbol{T}(-Ct)$. Therefore, we can calculate the velocity of the material point $s$ in the body frame as

$$\boldsymbol{u}_{S}(s,t)=C\boldsymbol{T}(-Ct)-C\boldsymbol{T}(s-Ct).$$

(31)

The material velocity with respect to the fluid at infinity is then

$$\boldsymbol{u}_{S}(s,t)+\boldsymbol{U}(t)+\boldsymbol{\Omega}(t)\times\boldsymbol{X}(s,t),$$

(32)

where $\boldsymbol{U}(t)$ and $\boldsymbol{\Omega}(t)$ are the translational and angular velocities of the flagellum.

Similarly to our treatment of an infinitely long flagellum in §2, we model the fluid velocity due to the flagellum’s motion as a linear superposition of velocities due to a force density (see equation (21)). Now the path of integration $\Gamma$ stands for the curve $\boldsymbol{X}(s,t)$ describing the instantaneous shape of the flagellum. The instantaneous swimming and angular velocities and the force density are calculated from the coupled integral equations for a no-slip boundary condition on the surface of the flagellum and the requirements of zero net force and torque on the flagellum,

$$\boldsymbol{u}_{S}(s,t)+\boldsymbol{U}(t)+\boldsymbol{\Omega}(t)\times\boldsymbol{X}=\frac{1}{4\pi\eta_{m}}\int_{\Gamma}\mbox{\boldmath$\alpha$}(\boldsymbol{X}-\boldsymbol{X}^{\prime})\bcdot\boldsymbol{f}({\boldsymbol{X}^{\prime}})\mathrm{d}{\boldsymbol{X}}^{\prime},$$

(33)

	$\displaystyle\int_{\Gamma}\boldsymbol{f}({\boldsymbol{X}})\mathrm{d}{\boldsymbol{X}}$	$\displaystyle=$	$\displaystyle 0,$		(34)
	$\displaystyle\int_{\Gamma}{\boldsymbol{X}}\times\boldsymbol{f}({\boldsymbol{X}})\mathrm{d}{\boldsymbol{X}}$	$\displaystyle=$	$\displaystyle 0,$		(35)

where ${\boldsymbol{X}}\equiv{\boldsymbol{X}}(s,t)$ and ${\boldsymbol{X}}^{\prime}\equiv{\boldsymbol{X}}(s^{\prime},t)$.

Similar to the approach discussed in §2, we solved the discretized version of these equations for the instantaneous velocities $\boldsymbol{\Omega}(t)$, $\boldsymbol{U}(t)$ and the force density $\boldsymbol{f}({\boldsymbol{X}})$. We normally calculated the angular and swimming velocities for about 60–80 snapshots per one period of oscillation and averaged them over one cycle of motion, $\langle\Omega\rangle=\frac{1}{T}\int_{0}^{T}\Omega(t)\mathrm{d}t\approx\frac{1}{N_{T}}\sum_{i=1}^{N_{T}}\Omega(t_{i})$ and $\langle\boldsymbol{U}\rangle=\frac{1}{T}\int_{0}^{T}\boldsymbol{R}(t)\cdot\boldsymbol{U}(t)\mathrm{d}t\approx\frac{1}{N_{T}}\sum_{i=1}^{N_{T}}\boldsymbol{R}(t_{i})\cdot\boldsymbol{U}(t_{i})$, where $\boldsymbol{R}(t)$ is the rotation operator that transforms the swimming velocity vector to the initial coordinate system $\boldsymbol{X}(s,0)$, which is motionless with respect to the fluid at infinity, and $N_{T}$ is the number of snapshots per one period.

During one cycle of motion, a flagellum of finite size moves in a tortuous fashion that involves pitching (rotation of the swimmer’s centerline with respect to the initial direction of wave propagation) and drifting (motion of the flagellum’s center of mass perpendicular to the net swimming direction), see figure 4.

Our computations predict that $\langle\Omega\rangle=0$, and a flagellum of finite size swims in a straight line at an angle $\gamma$ with respect to $(-x)$-axis, where $\gamma=\arctan(\langle U_{y}\rangle/|\langle U_{x}\rangle|)$. The magnitude and sign of $\gamma$ depend on the flagellum’s contour length, the wave amplitude and the phase constant $\phi_{0}$. For particular values of $\phi_{0}$, a flagellum of a given contour length can assume an even or odd configuration at time $t=0$. In the even configuration (as illustrated in figure 5b) the flagellum has reflection symmetry with respect to the vertical dashed line that passes through the flagellum’s center of mass. In the odd configuration the shape of the flagellum has point symmetry about the center of the flagellum (marked by cross hairs in figure 5a).

Koehler et al. (2012), in work on the swimming of finite-length flagella in a Newtonian 3D fluid, used symmetry arguments to prove that a flagellum that starts its motion from an even configuration swims along its centerline in the direction opposite to initial wave propagation, and does not drift away from the $(-x)$-direction. Koehler et al. (2012) pointed out that the mirror reflection about the vertical line is equivalent to the time reversal, with the time-reversed swimming velocity $\boldsymbol{U}_{-t}=-\boldsymbol{U}_{t}$. Therefore, the instantaneous swimming velocity in an even configuration must be identical to the mirror image of the time-reversed velocity. This condition requires the $y$- component of the swimming velocity be equal to zero. Thus, a flagellum starting its motion from an even configuration, has vanishing mean drift over one cycle of motion. A flagellum starting its motion in an odd configuration, will reach an even configuration after a quarter of a period. The reorientation angle acquired by the flagellum by $t=T/4$ determines the flagellum’s swimming direction. The reorientation angle for a flagellum starting its motion in an odd configuration is the largest since it takes the longest amount of time for such a flagellum to reach an even configuration. In figure 5 we plot $\gamma$ as a function of the flagellum’s contour length $N_{\Lambda}\equiv L/\Lambda$ for various values of the phase constant $\phi_{0}$. A similar swimming pattern was reported by Peng et al. (2016) in the work on flagella locomotion in granular media.

In figure 6a we plot the swimming speed averaged over one period of oscillations, $U=|\langle\boldsymbol{U}\rangle|$, as a function of a dimensionless parameter $bq$ for the flagellum contour length equal to one arcwise wavelength, $N_{\Lambda}=1$. Two competing mechanisms influence the swimming speed of the flagellum. On the one hand, larger values of $bq$ correspond to steeper angles between the flagellum and the direction of wave propagation and, therefore, a stronger propulsion force. On the other hand, for larger $bq$ values the segments of the flagellum come closer to each other. The hydrodynamic interactions between the segments tend to slow down the swimmer. The hydrodynamic interactions are stronger in the limiting case of 2D hydrodynamics (small $\lambda/\ell_{S}$ ratios) due to a slow, logarithmic spatial decay rate of the flow field. In the opposite limit of large $\lambda/\ell_{S}$ ratios, our calculations do not reproduce Higdon’s results for the swimming speed in a purely 3D fluid (Higdon, 1979). Being qualitatively similar to Higdon’s prediction, our calculations show much larger swimming speeds for $\lambda/\ell_{S}\gg 1$. As we discussed at the end of §2, the incompressibility of the membrane sets a constraint on the fluid dynamics that leads to an effective drag anisotropy that grows logarithmically as a function of $\lambda/\ell_{S}$ for $\lambda/\ell_{S}\gg 1$. The enhanced drag anisotropy in a quasi-2D membrane is responsible for larger swimming speeds in quasi-2D membranes (in comparison with pure 2D or 3D fluids).

In figure 6b we plot the swimming speed as a function of the scaled flagellum length $N_{\Lambda}$ for $bq=1$. For $N_{\Lambda}<1$ the flagellum performs large yawing motion that is inefficient for swimming (see figure 4a). For larger values of $N_{\Lambda}$ the long-range hydrodynamic interactions taper off the growth of the swimming speed, and the speed approaches the values found for an infinitely long flagellum (shown as dotted black horizontal lines in figure 6b). The ‘bumps’ in the curves reflect smaller yawing of the flagellum for some values of $N_{\Lambda}$.

To find the flagellum motion that is optimal in terms of the power consumption, we calculated the swimming efficiency. As discussed in (Koehler et al., 2012), there are multiple efficiency metrics. Here we calculated the efficiency as the ratio of the power required to pull the flagellum through the fluid at its average swimming speed, $F_{\rm pull}U=(1/\mu_{\|})U^{2}$, to the average power $\langle P\rangle$ consumed by the swimmer over one period of motion,

$$\eta=\frac{(1/\mu_{\|})U^{2}}{\langle P\rangle}.$$

(36)

Here $\mu_{\|}$ is the flagellum mobility for the translational motion along the direction of wave propagation, averaged over one period. The average power consumption is

$${\langle P\rangle}=\frac{1}{T}\int_{0}^{T}\mathrm{d}t\int_{\Gamma}\mathrm{d}{\boldsymbol{X}}\boldsymbol{f}({\boldsymbol{X}},t)\cdot\boldsymbol{u}_{S}({\boldsymbol{X}},t).$$

(37)

We calculated the mobility and the power consumption numerically using the BEM.

Our calculated flagellum efficiency is in qualitative agreement with Higdon (1979). In figure 7a we plot the efficiency as a function of the amplitude $bq$ for a few values of the scaled wavelength $\lambda/\ell_{S}$ for a flagellum of length $N_{\Lambda}=1$. For smaller values of $bq$ the segments of the flagellum have small angles with respect to the direction of wave propagation and produce a weak thrust. For larger values of $bq$ the flagellum ‘shrinks’ along the $x$-axis, and the stronger interference between the segments of the flagellum leads to a decrease in efficiency. The efficiency increases with $\lambda/\ell_{S}$ due to the reduced role of long-range hydrodynamics on length scales exceeding the Saffman length $\ell_{S}$. The maximum efficiency falls at $bq=$ 1.4 – 1.8.

In figure 7b we plot the efficiency as a function of the flagellum length $N_{\Lambda}$ for a wavelength $\lambda/\ell_{S}=1$. For small values of $N_{\Lambda}$ the swimming of the flagellum is inefficient due to an excessive yawing motion and a weak overall thrust (see figure 4a). The efficiency reaches a maximum at $N_{\Lambda}$ = 0.8 – 1.4 and then decreases with further growth of $N_{\Lambda}$ due to interference between the crests of the flagellum. The interference is stronger for larger amplitudes $bq$ since the crests are closer to each other, and the efficiency drops off more abruptly from its optimal value for larger values of $bq$. The secondary maxima correspond to a smaller yawing and larger propulsion for various values of $N_{\Lambda}$. Figure 4 demonstrates that the flagellum travels a considerable distance along the $y$-axis while making moderate overall progress along the $x$-axis.

In the following section we discuss an alternative computational approach to finding the swimming velocities using the Lorentz reciprocal theorem.

Finding the analytical solution for the swimming velocity can be a daunting task. One of the major difficulties is that one needs to solve Stokes equations with time dependent no-slip boundary conditions on the surface of the swimmer. Stone & Samuel (1996) offered an elegant way to find the swimming velocity using the Lorentz reciprocal theorem (LRT) (Happel & Brenner, 1965). For a fluid in 3D the Lorentz reciprocal theorem states that if there are two solutions to Stokes equations and the incompressibility condition with the velocity fields and the stress tensors $(\boldsymbol{v},\ \boldsymbol{\sigma})$ and $(\boldsymbol{v}^{\prime},\ \boldsymbol{\sigma}^{\prime})$, respectively, that satisfy the same boundary conditions at infinity, than for a volume of fluid $V$ bounded by surface $S$, we have

$$\oint_{S}\boldsymbol{v}\cdot\boldsymbol{\sigma}^{\prime}\cdot\boldsymbol{n}\,\mathrm{d}S=\oint_{S}\boldsymbol{v}^{\prime}\cdot\boldsymbol{\sigma}\cdot\boldsymbol{n}\,\mathrm{d}S,$$

(38)

where $\boldsymbol{n}$ is the outward normal to the surface $S$.

Here we formulate the LRT approach for a finite swimmer confined to a quasi-2D membrane. Let $\boldsymbol{v}$ and $\boldsymbol{\sigma}$ be the membrane velocity and the stress fields for the swimming problem. These velocity and stress fields are solutions of the Stokes equations and the incompressibility condition, equation (4). They also satisfy the conditions of zero net force and torque on the swimming body, equations (7) and (8). For the reciprocal solution of equation (4) we choose the membrane velocity ${\boldsymbol{v}}^{\prime}$ and stress ${\boldsymbol{\sigma}}^{\prime}$ fields due to an inactive object of the same shape as the swimmer and being dragged as a solid body with constant translational velocity ${\boldsymbol{U}}^{\prime}$.

When the condition $\bnabla\cdot\boldsymbol{\sigma}=0$ is relaxed (see equation (6)), a more general form of the Lorentz reciprocal theorem is (Kim & Karrila, 1991)

$\displaystyle\oint_{S}{\boldsymbol{v}}^{\prime}\cdot(\boldsymbol{\sigma}\cdot\boldsymbol{n})\mathrm{d}S-\int_{V}{\boldsymbol{v}}^{\prime}\cdot(\bnabla\cdot\boldsymbol{\sigma})\mathrm{d}V=\oint_{S}\boldsymbol{v}\cdot({\boldsymbol{\sigma}}^{\prime}\cdot\boldsymbol{n})\mathrm{d}S-\int_{V}\boldsymbol{v}\cdot(\bnabla\cdot{\boldsymbol{\sigma}}^{\prime})\mathrm{d}V,$

(39)

where $V$ is the swimmer’s volume bounded by the surface $S$. Here we treat the volume occupied by the swimmer as being equivalent to the fluid domain of the same shape as the swimmer’s and having the same velocity distribution as that of the material particles of the swimmer.

Let us consider the first term on the LHS of equation (39):

$\displaystyle\oint_{S}{\boldsymbol{v}}^{\prime}\cdot(\boldsymbol{\sigma}\cdot\boldsymbol{n})\mathrm{d}S={\boldsymbol{U}}^{\prime}\cdot\int_{S_{w}}\mathrm{d}\boldsymbol{F},$

(40)

where we took into account that ${\boldsymbol{v}}^{\prime}={\boldsymbol{U}}^{\prime}$ is a constant vector at the surface of the domain (uniform translation). Also, since $\boldsymbol{\sigma}$ does not have $z$- components, only $\boldsymbol{\sigma}\cdot\boldsymbol{n}\mathrm{d}S=\mathrm{d}\boldsymbol{F}$ on the curvy wall of the domain, $S_{w}$ (see figure 8), will make a non-zero contribution.

Taking into account equation (6), the second term on the LHS of equation (39) can be rearranged as

$\displaystyle-\int_{V}{\boldsymbol{v}}^{\prime}\cdot(\bnabla\cdot\boldsymbol{\sigma})\mathrm{d}V={\boldsymbol{U}}^{\prime}\cdot\int_{V}\frac{2\boldsymbol{f}}{h}\mathrm{d}V={\boldsymbol{U}}^{\prime}\cdot\int_{S_{t,b}}2\boldsymbol{f}\,\mathrm{d}S,$

(41)

where we take into account that the traction forces $2\boldsymbol{f}$ due to the fluid flows in the surrounding fluid act on the flat sides of the domain $S_{t,b}$ (see equation (6) and figure 8), and $dV=hdS$. Therefore, the LHS of equation (39) becomes:

$\displaystyle{\boldsymbol{U}}^{\prime}\cdot\left(\int_{S_{w}}\mathrm{d}\boldsymbol{F}+\int_{S_{t,b}}2\boldsymbol{f}\,\mathrm{d}S\right)={\boldsymbol{U}}^{\prime}\cdot{\boldsymbol{F}}=0,$

(42)

where ${\boldsymbol{F}}$ is the net force on the swimmer and is equal to zero.

Similarly, for the terms on the RHS of equation (39) we have

	$\displaystyle\oint_{S}\boldsymbol{v}\cdot({\boldsymbol{\sigma}}^{\prime}\cdot\boldsymbol{n})\mathrm{d}S-\int_{V}\boldsymbol{v}\cdot(\bnabla\cdot{\boldsymbol{\sigma}}^{\prime})\mathrm{d}V$	$\displaystyle=$	$\displaystyle\int_{S_{w}}\boldsymbol{v}\cdot({\boldsymbol{\sigma}}^{\prime}\cdot\boldsymbol{n})\mathrm{d}S+\int_{S_{t,b}}\boldsymbol{v}\cdot(2{\boldsymbol{f}}^{\prime})\,\mathrm{d}S$		(43)
		$\displaystyle=$	$\displaystyle\oint_{S}\boldsymbol{v}\cdot\mathrm{d}{\boldsymbol{F}}^{\prime},$

where in the last line of equation (43) we merged two terms into one integral over the total surface of the swimmer, and $\mathrm{d}{\boldsymbol{F}}^{\prime}$ denotes an elementary traction force on the inactive ‘swimmer’ being dragged with constant velocity. Thus, the Lorentz reciprocal relation, equation (39), assumes a compact form,

$$0=\oint_{S}\boldsymbol{v}\cdot\mathrm{d}{\boldsymbol{F}}^{\prime}.$$

(44)

Decomposing the surface velocity of the swimmer into the translational $\boldsymbol{U}(t)$ and the surface disturbance $\boldsymbol{u}_{S}(t)$ velocities, $\boldsymbol{v}(t)=\boldsymbol{U}(t)+\boldsymbol{u}_{S}(t)$, we rewrite the Lorentz reciprocal relation in the form

$${\boldsymbol{F}}^{\prime}(t)\cdot\boldsymbol{U}(t)=-\oint_{S(t)}\boldsymbol{u}_{S}\cdot\mathrm{d}{\boldsymbol{F}}^{\prime},$$

(45)

similar to the equation derived by Stone & Samuel (1996) for a swimmer in a 3D fluid. In equation (45) the integration is performed over the instantaneous surface area of the swimmer. The generalization of (45) for the motion that involves rotation is

$${\boldsymbol{F}}^{\prime}(t)\cdot\boldsymbol{U}(t)+{\boldsymbol{L}}^{\prime}(t)\cdot\boldsymbol{\Omega}(t)=-\oint_{S(t)}\boldsymbol{u}_{S}\cdot\mathrm{d}{\boldsymbol{F}}^{\prime},$$

(46)

where ${\boldsymbol{L}}^{\prime}(t)$ is the torque applied to the inactive inclusion and $\boldsymbol{\Omega}(t)$ is the swimmer’s angular velocity.

The Lorentz reciprocal relation, equation (46), is particularly useful for computation of the swimmer’s translational and rotational velocities, $\boldsymbol{U}(t)$ and $\boldsymbol{\Omega}(t)$, when the stress tensor of the reciprocal problem (motion of an inactive body) is known. Unfortunately, it is also difficult to solve the reciprocal problem analytically for an inclusion of an arbitrary shape in a quasi-2D membrane, since the coupling with the bulk fluid makes the problem essentially three-dimensional. When the reciprocal solution is not available, equation (46) can serve as an alternative computational path for finding the swimming velocity.

As a test of equation (46), we found the swimming velocities of a finite inextensible flagellum as described in §3. We solved the reciprocal problem numerically by finding the force densities $\boldsymbol{f}({\boldsymbol{x}},t)$ for the uniform rotation of the flagellum about the $z$-axis and for the translational motion of the flagellum along the $x$- and $y$- axes for multiple flagellum conformations corresponding to various time instants of the swimming cycle. We then solved the resulting system of three equations (46) for the instantaneous angular velocity and $x$- and $y$- components of the translational velocity and found the average swimmer’s speed over one period of oscillations. In figure 6b circles superimposed on the curves show the swimming velocities obtained using the Lorentz equation (46).

We have studied analytically and computationally the locomotion of microscopic organisms confined to a plane of a thin fluid membrane embedded in a bulk fluid of different viscosity. In our model the membrane is sufficiently thin, with material particles moving only in the plane of the membrane (the motion in the perpendicular direction is forbidden).

The presence of the bulk fluid allows the introduction of a hydrodynamic length scale, the Saffman length, that controls the energy exchange between the membrane and the surrounding fluid. By varying the Saffman length, we make our model continuously vary between a pure 2D system (large Saffman length) and a quasi-2D system (small Saffman length). The hydrodynamic flows in the quasi-2D membrane have features of both 3D and 2D hydrodynamics. We show that a flagellated swimmer in a viscous film (Saffman length smaller than swimmer characteristic length scale) swims faster than the same swimmer in a 3D fluid. The speedup comes from the effectively larger perpendicular drag coefficient, which arises from the incompressibility of the membrane. On the other hand, a circular squirmer, whose propulsion mechanism does not employ the local drag anisotropy, slows down for smaller Saffman lengths (in comparison with the squirmer’s radius).

The coupling of the membrane with the bulk fluid makes the problem three-dimensional and quite difficult for analytical treatment. We developed numerical schemes based on the boundary element method and the Lorentz reciprocal theorem. We show how the Lorentz reciprocal theorem can be used to simplify the computation of swimming speed, especially for swimmers such as the squirmer that do not change shape during a stroke.

While we considered the minimal models of a flagellated swimmer and of a squirmer, our approach can be generalized to other swimmers’ geometries and swimming stokes.

Acknowledgements We thank M. Moelter for helpful discussions. This work has been supported by a Cottrell College Science Award (TK) and National Science Foundation Grant No. 1437195 (TRP). (CA) gratefully acknowledges support from a Frost Undergraduate Student Research Award.