Analytical Study for Swimmers in a Channel

(1)

vol. 881, pp. 365–383. c Cambridge University Press 2019 doi:10.1017/jfm.2019.751

365

Analytical study for swimmers in a channel

A. Farutin^1,†, H. Wu¹, W.-F. Hu², S. Rafaï¹, P. Peyla¹, M.-C. Lai³ and C. Misbah¹

1Univ. Grenoble Alpes, CNRS, LIPhy, F-38000 Grenoble, France

2Department of Applied Mathematics, National Chung Hsing University, 145, Xingda Road, Taichung City 402, Taiwan

3Department of Applied Mathematics, National Chiao Tung University, 1001, Ta Hsueh Road, Hsinchu 300, Taiwan

(Received 29 March 2019; revised 25 June 2019; accepted 25 August 2019)

There is an overabundance of microswimmers in nature, including bacteria, algae, mammalian cells and so on. They use flagellum, cilia or global shape changes (amoeboid motion) to move forward. In the presence of confining channels, these swimmers exhibit often non-trivial behaviours, such as accumulation at the wall, navigation and so on, and their swimming speed may be strongly influenced by the geometric confinement. Several numerical studies have reported that the presence of walls either enhances or reduces the swimming speed depending on the nature of the swimmer, and also on the confinement. The purpose of this paper is to provide an analytical explanation of several previously obtained numerical results. We treat the case of amoeboid swimmers and the case of squirmers having either a tangential (the classical situation) or normal velocity prescribed at the swimmer surface (pumper).

For amoeboid motion we consider a quasi-circular swimmer which allows us to tackle the problem analytically and to extract the equations of the motion of the swimmer, with several explicit analytical or semi-analytical solutions. It is found that the deformation of the amoeboid swimmer as well as a high enough order effect due to confinement are necessary in order to account for previous numerical results.

The analytical theory accounts for several features obtained numerically also for non-deformable swimmers.

Key words: swimming/flying, micro-organism dynamics

1. Introduction

Microorganisms swimming under low Reynolds number conditions use different locomotion strategies. Many studied microswimmers use cilia or flagella (Drescher et al. 2009; Lauga & Powers 2009; Guasto, Johnson & Gollub 2010; Garcia et al.

2011; Saintillan & Shelley 2012; Kantsler et al. 2013). Prototypical examples are Chlamydomonas reinhardtii, Escherichia coli and Paramecium. Another mode of locomotion, prevalent among eukaryotic cells, albeit much less studied theoretically, is amoeboid motion. This mode is characterized by large deformations of the cell.

† Email address for correspondence: alexandr.farutin@univ-grenoble-alpes.fr

https://www.cambridge.org/core. National Chiao Tung University, on 28 Oct 2019 at 01:54:34, subject to the Cambridge Core terms of use, available at https://www.cambridge.org/core/terms.

(2)

A prototypical example is Eutreptiella gymnastica (Throndsen 1969). The name amoeboid is generically used for any motion based on body deformation (be it crawling on a substrate or swimming). Several studies have reported that amoeba such as Dictyostelium (Bae & Bodenschatz 2010; Barry & Bretscher 2010), but also leucocytes (Bae & Bodenschatz 2010; Aoun et al. 2019) or even cancer cells (Pinner

& Sahai 2008) can swim. Thus, the swimming of these types of microorganisms without adhesion assistance is now gaining more and more interest in the biological literature (Hawkins et al. 2009; Bergert et al. 2015; Liu et al. 2015). Understanding amoeboid swimming in a simple fluid has recently incited several theoretical studies (Shapere & Wilczek 1987; Avron, Gat & Kenneth 2004; Ohta & Ohkuma 2009;

Alouges, Desimone & Heltai 2011; Hiraiwa, Shitara & Ohta 2011; Vilfan 2012;

Farutin et al. 2013; Loheac, Scheid & Tucsnak 2013; Wu et al. 2015, 2016).

The literature on the wall effect on swimming is quite abundant. The presence of a wall may add several features to the swimming of microorganisms: it may lead to (i) circular trajectories (Smith et al. 2009), (ii) accumulation of swimmers at the wall (Ezhilan & Saintillan 2015), (iii) the existence of a stable position of the swimmer at a given distance from the wall (Giacché, Ishikawa & Yamaguchi 2010; Shum, Gaffney & Smith 2010), (iv) rich dynamics of the swimmer orientation and scattering (Crowdy & Or 2010; Lushi, Kantsler & Goldstein 2017) and so on. A classical way of representing the presence of the walls is via multipolar representation. The first natural step to study swimming in a confined geometry is to represent the swimmer by leading-order terms in far field expansions, which is the Stokeslet dipole. It has been shown (Spagnolie & Lauga 2012) that the presence of the wall may strongly influence the swimming behaviour.

Our study is focused here on the determination of the swimming speed between two walls, where several numerical observations have not yet been elucidated. Several groups have studied the problem of the influence of confining walls on swimming.

Felderhof (2010) has reported analytically on a Taylor-like swimmer, showing that the speed is enhanced with confinement. Zhu, Lauga & Brandt (2013) have considered numerically the case of a squirmer to show that (when only tangential surface flow is included) the velocity decreases with confinement. When considering normal deformation, they found an increase of velocity with confinement. Liu et al. (2015) have studied numerically a helical flagellum in a tube and found that, except for a small range of tube radii, the swimming speed, when the helix rotation rate is fixed, increases monotonically as the confinement becomes tighter. Acemoglu & Yesilyurt (2014) have adopted a similar model but, besides the flagellum, their swimmer possesses a head and they found a decrease of velocity with confinement. Bilbao et al. (2013) have dealt numerically with a model inspired by nematode locomotion and found that it moves faster due to walls. Ledesma-Aguilar & Yeomans (2013) have analysed a dipolar swimmer in a rigid or elastic tube and found a speed enhancement due to the walls. Recently, we have studied numerically the case of amoeboid swimmers (Wu et al. 2015, 2016). The velocity was found to decrease slightly with increasing confinement for very wide channels, then to significantly increase for intermediate confinement and finally to decrease to zero as the channel width was further decreased. Understanding of this complex behaviour of amoeboid swimmers is our first focus. It will be shown here that the deformation as well as higher-order effects due to confinement play a decisive role here. We are not aware of any previous study taking these effects into account. We will also consider the case of a squirmer and pumper (undeformable swimmer on the surface of which the normal velocity is prescribed – the swimmer may be thought of as a porous medium pumping fluid) studied numerically (Zhu et al. 2013).

(3)

We build a systematic analytical framework in order to derive the expression of the swimming speed of a microswimmer under confinement. We consider the case of circular swimmers with fixed shape (with tangential or normal flow) and the case of amoeboid swimmers. The obtained analytical expressions agree quantitatively with numerical results for weakly and moderately confined swimmers. For amoeboid swimmers, we observe that the leading effect of confinement (inversely proportional to the square of the channel width) is anomalously small due to an almost-perfect counterbalance of two effects: (i) a wall-induced decrease of the swimming velocity due to a source-dipole term in the multipole expansion of the flow generated by the swimmer and (ii) an increase of the swimming velocity due to swimming in an effective extensional flow which is induced by the walls in reaction to the stresslet of the swimmer. The next-order effect (inversely proportional to the fourth degree of the channel width) has a positive coefficient and is responsible for the increase of the swimming velocity at moderate confinements. For the squirmer case, we find that confinement decreases the speed when tangential flow is considered (the classical case). In the case of normal flow, the speed is enhanced. These results agree with numerical results obtained in three dimensions (Zhu et al.2013), despite the fact that our analytical work is developed in two dimensions. This highlights the fact that the dimensionality is irrelevant.

In the case of amoeboid swimmer the treatment is more complex than the squirmer case, due to the free-boundary character. For the sake of analytical tractability, we consider that the swimmer has a shape which is not far from a circle. This allows us to perform a systematic expansion in power series of a small parameter, namely the deviation from a circle. The calculations are simplified by using complex notations for vector variables and representation by analytical functions for multipole expansions.

2. The model

We first consider amoeboid motion as a general case. The case of undeformable swimmers will follow from the general consideration as a particular limit. We use the two-dimensional model introduced in the numerical studies (Wu et al. 2015, 2016) in order to facilitate direct comparison of the results. The system set-up is drawn schematically in figure 1. The swimmer is modelled as an inextensible membrane enclosing a viscous fluid. A time-dependent distribution of active forces normal to the membrane is applied by the swimmer to the inner and outer fluids. The active forces are chosen to be periodic in time, where the period is defined as one stroke cycle.

The membrane deforms due to the active forces, which can lead to a net displacement of the swimmer after one stroke, provided the deformation is not time reversible. A translational swimming motion can be achieved by repeating the strokes. We limit this work to the case of swimming along the axis of the channel. It was observed in studies (Wu et al. 2015, 2016) that such centred swimming is not stable and the actual stable mode involves navigation from one wall to the other and back in addition to translation along the channel axis. Nevertheless, the velocity of the translation along the channel axis was found to be close to that in the centred mode.

The size of the swimmer is characterized by the area of the enclosed fluid A0, which is conserved due to the impermeability of the membrane and the incompressibility of the fluid inside the swimmer. Another geometric parameter is the perimeter of the swimmer L0, which we assume to be constant due to the inextensibility of the membrane enclosing the swimmer. These two parameters define a non-dimensional number τ ≡ 4πA0/L²0, which is called the reduced area of the swimmer. A related

(4)

-W/2 W/2 y

0 x

FIGURE 1. (Colour online) Amoeboid swimmer in a channel.

non-dimensional number is the excess perimeter (counted from the circular shape) defined as

Γ ≡ L₀ 2π√

A0/π−1 =τ⁻¹^/2−1. (2.1) We use the excess perimeter Γ as a small parameter in which the swimming velocity is expanded. The case Γ = 0 corresponds to a perfectly circular swimmer, while positive values of Γ specify the extent of the deformation amplitude.

The swimmer is suspended in a Newtonian fluid of viscosity η. We assume the interior region of the swimmer to be filled with a fluid having the same viscosity in order to simplify the derivation. Nevertheless, the viscosity of the inner fluid does not affect the swimming velocity in the fully saturated regime which we treat here (i.e.

when the active force is so large that the swimmer has enough time to fully develop deformation amplitude available to it before starting a new elementary stroke), as discussed below. The fluids inside and outside the swimmer obey the Stokes equations

−∇P(r) + η∇²u(r) = 0,

∇ ·u(r) = 0,

(2.2) where r is the position vector, P is the pressure and u is the velocity.

The swimmer is confined between two rigid walls separated by a distance W (figure 1). We impose the no-slip boundary condition at the walls. The confinement parameter of the problem is defined as Cn = 2R/W, where R = √

A0/π is the characteristic size of the swimmer. We use Cn as the second small parameter in which the swimming velocity is expanded. The case C_n=0 corresponds formally to an unconfined swimmer. We place the swimmer in the centre of the channel with swimming velocity parallel to the walls. The wall direction is used as the x coordinate axis and the swimmer position as the origin (as shown in figure 1).

The motion of the swimmer is actuated by forces applied from the membrane to the fluid. The total force density (force per unit area) at the membrane is composed of an active part and a passive part, and can be written as (Wu et al. 2016)

f = f_a−ζcn +∂ζ

∂st + f₀+ftt, (2.3) where f_a is the density of the active force, which depends on position and time. The exact expression of the density of the active force is chosen below; ζ is a Lagrange

(5)

multiplier that enforces local inextensibility of the membrane, c is the curvature, t is the unit tangent vector and s is the arc length. The constants f0 and f_t are added to satisfy the condition of the total force and the total torque exerted by the membrane being equal to zero,

I

f ds = 0, I

r × f ds = 0. (2.4a,b)

The symmetry of the problem considered here dictates f_0y=ft=0.

We use the same expression of the active force as in Wu et al. (2015, 2016), f_a(s, t) = 2A[cos(ωt) cos(2α) − sin(ωt) cos(3α)]n, (2.5) where t is time, ω is the stroke frequency, α ≡ 2πs/L0 is the rescaled arc length on the swimmer membrane and A controls the amplitude of the active stresses.

Local incompressibility of the membrane of the swimmer implies that α serves as a reference coordinate of material points on the membrane of the swimmer. Indeed, as the shape of the swimmer evolves with time, the distance along the membrane between any two material points remains constant.

The non-dimensional saturation number S ≡ A/(ηω) expresses the amplitude of the active stresses compared to the viscous ones. The saturation number S can be viewed as a ratio of two time scales: the characteristic time of shape response under the action of active stresses η/A and the characteristic time over which the active forces change 1/ω. Large values of S correspond to saturated regime, in which the shape adapt almost instantaneously to the changes of the active stresses. In this study we consider a formal limit of S → ∞, which closely resembles the quasi-saturated case of S = 10 used to obtain most of the results in Wu et al. (2015, 2016). An estimate for S was given in Ranganathan, Farutin & Misbah (2018) for mammalian cells and was found to be of order one. This is a very rough estimate since the cortex viscosity (which is the most relevant one; see Ranganathan et al. (2018)) is poorly documented.

In our previous numerical study (Wu et al. 2016) we found that the swimming speed is close to saturation even at S = 1. In view of this we have opted for the limit where S is large enough, a limit which lends itself to a relatively simpler analytical tractability.

The solution to the problem is obtained as a truncated expansion in powers of Γ¹^/2 (we show below the legitimacy of this choice) and Cn. The convergence of the expansions for given values of Γ and Cⁿ is verified by comparison to full numerical simulations of swimming in confined geometry. The value S = 5000 was used in numerical simulations.

3. Analytical technique

3.1. Complex notation

This section provides the main steps of the analytical solution of the problem. We present explicitly the solution to the leading order in powers of Γ and Cn. Higher- order expansions are obtained by following the same procedure but we only provide the final expression for the time-averaged swimming velocity because the intermediate expressions are too long.

It is convenient to view two-dimensional vectors as complex numbers, which we denote by a hat. For example, ˆr = rx+iry. All vector fields on the swimmer membrane can then be represented by periodic complex-valued functions of α. These periodic

(6)

functions are conveniently represented by Fourier series. For example, the shape of the swimmer is parametrized as

rˆ(α) =

∞

X

k=−∞

rke^ik^α. (3.1)

3.2. Quasi-static approximation

The first step is to apply the quasi-static approximation, corresponding to the limit S → ∞. The stresses in fluids adjacent to the membrane of the swimmer are balanced by the force density acting from the membrane. The velocity and pressure gradients inside the fluids are controlled by the deformation rate of the membrane, which, in turn, is controlled by the stroke frequency in the limit of large S. Therefore, the stresses in fluids are small compared to A in the limit of large S, except for the constant pressure difference 1P between the inner and outer fluids. Taking the limit S → ∞yields for the force balance at the membrane ( f + 1Pn)/A = 0, where n is the normal to the membrane. In other words, the non-hydrostatic stress is irrelevant in the limit of large active forces. This can also be formally seen by appropriate rescaling of variables. Indeed, rescaling time by ω yields in the stress balance a factor ηω for the amplitude of the non-hydrostatic stress, to be compared to A for the active force. The ratio between the latter and the former provides S (which is taken to be large enough in our study). The normal and tangential projections of the simplified force balance equation read, respectively,

fa−ζc + f0nx+1P = 0, (3.2)

∂ζ

∂s +f0tx=0, (3.3)

where we have omitted the denominator, which is irrelevant for further derivation, in order to simplify the notation. Since tx=drx/ds, equation (3.3) can be integrated as

ζ + f⁰r_x=ζ⁰, (3.4)

where ζ0 corresponds to the homogeneous part of the tension.

3.3. Small-deformation approximation

We parametrize the shape of the swimmer by the tangent direction as a function of α because this allows us to ensure the local inextensibility of the membrane automatically,

dˆr dα =iL0

2πeⁱ^αeⁱ^φ(α), (3.5)

where φ(α) is a real, time-dependent function used to parametrize the shape of the swimmer. Indeed, the arc element on the swimmer contour is written as ds = |dr| = (|dr|/dα) dα. Consequently, the absolute value of dˆr/dα, which is equal to L0/(2π) in (3.5), specifies the local stretching of the swimmer boundary with respect to a circle of radius 1. This shows that parametrization (3.5) allows us to follow not just the shape of the swimmer but also the position of material points on the boundary of a swimmer with a locally inextensible membrane. The argument of dˆr/dα, which is equal to π/2 + α + φ(α), specifies the angle which the tangent vector of the swimmer boundary

(7)

makes with the x axis. This is because a unit vector with coordinates(cos α, sin α) is conveniently written as eⁱ^α in the complex notation used here. We have split off the eⁱ^α exponent in (3.5) in order for (3.5) to represent a circle of perimeter L0 for φ = 0. The prefactor i could actually be adsorbed in the exponent but we write it here explicitly in order to simplify the final expressions: indeed, imposing the y → −y symmetry of the problem is equivalent to requiringφ to be an odd function of α. Consequently, φ can be expanded as

φ(α) =

∞

X

k=−∞

φksin(kα). (3.6)

The small-deformation approximation relies on the amplitudes φ^k being small, which allows us to make formal expansions of all shape-related equations. We introduce an expansion parameter which expresses the smallness of each particular term. We show below that = O(Γ¹^/2). Because the active force (2.5) contains only the second and the third harmonics of α, we assume φ² and φ³ to be of order O(). This ansatz is validated by the consistency of the obtained expansions. The amplitudesφ^k for k from 4 to 6 are of order O(²), and so on. The amplitude φ1 is a special case because it can be expressed through the other amplitudes by imposing the consistency condition

Z 2π 0

dˆr

dαdα = 0. (3.7)

Expanding (3.7) in powers of φ yields iL0

2π Z 2π

0

eⁱ^α

1 + iφ(α) −1 2φ(α)²

dα + O(³) =iL0

2π

−πφ1−πφ²φ³ 2

+O(³) = 0.

(3.8) Equation (3.8) shows that φ¹ = −φ²φ³/2 + O(³). The amplitudes r^k in (3.1) are expressed through amplitudes φk by integrating equation (3.5). The O() terms are written as

r−2= L0

8πφ3+O(²), r−1= L0

4πφ2+O(²), r3= L0

12πφ2+O(²), r4= L0

16πφ3+O(²).











(3.9)

The other coefficients rk are of higher order. The amplitude r1 has to be expanded to the next order in order to obtain the leading order of the area inside the swimmer

r1= L0

2π

1 −φ2²+φ3²

4

+O(³). (3.10)

The constant part r0 of expansion (3.1) cannot be deduced from (3.5) because the shape is translationally invariant. For simplicity, we set r₀=0, effectively choosing the origin to be located in the centre of mass of the membrane of the swimmer. The actual displacement is calculated below by time integration of the velocity of the centre of mass of the membrane. This, together with (3.4), allows us to express the Lagrange multiplier ζ through φk, ζ0 and f0

ζ = ζ0−L0f0cosα

2π +O(²), (3.11)

where we assume that f0=O(), as is validated below.

(8)

The expression for ζ is then substituted in (3.2), where the curvature c can be conveniently written as

c =2π L0

1 +dφ

dα

. (3.12)

The consistency of the equations requires us to assume at this point that ζ⁰ scales as O(⁻¹). Such a non-analytical behaviour of the membrane tension is also encountered for quasi-spherical inextensible membranes in three dimensions (Danker et al. 2007).

The amplitudes φ^k are obtained by solving the kth Fourier harmonic (the coefficient of e^ik^α) in (3.2),

φ²=Acos(ωt) 4πζ0

+O(²), φ³= −Asin(ωt) 6πζ0

+O(²). (3.13a,b) The first harmonic of the same equation yields the expansion of f0,

f0=L0cos(t) sin(ωt) 48πζ0

. (3.14)

Had we supposed that the tensionζ0 scales differently, we would then have obtained that the leading shape amplitudes φ2 and φ3 vanish, which is absurd (no symmetry constraint imposes vanishing of φ² and φ³). The zeroth harmonic can be used to calculate the pressure difference 1P but this quantity is not needed for the subsequent calculations. At this point the only quantity remaining unknown is the value ofζ0. The expression of ζ⁰ is obtained by calculating the area inside the swimmer as

A0 = 1 2Im

Z 2π 0

rˆdˆr^∗

dα dα = L₀² 4π

1 − 9

16φ²2−2

3φ3²+O(³)

= L²₀ 4π

1 − L²₀ 64π²ζ0²

sin²(ωt) − L²₀ 24π²ζ0²

cos²(ωt) + O(³)

. (3.15)

Substituting equation (3.15) into (2.1) yields the expansion of ζ0 in powers of Γ with time-dependent coefficients.

Γ =1 3φ2²+ 9

32φ3²+O(³) =L²₀[5 cos²(ωt) + 3]

384π²ζ0²

+O(³). (3.16) This equation clearly shows that the deviation from the circular shape, measured by φk, is indeed of order Γ¹^/2, as anticipated. Equation (3.16) has two solutions for ζ0

as a function of Γ , only one of which is stable against shape perturbation for finite S.

This can be seen by referring to the shape evolution equation (4) (Farutin et al. 2013) (valid for small Γ ) where linear stability analysis shows that the steady-state solution is unstable for ζ⁰< 0. A physical argument is that ζ⁰< 0 corresponds to the situation where the membrane is under compression, which leads to buckling.

The solution reads

ζ0=L₀[18 + 30 cos²(ωt)]¹^/2

48πΓ¹^/2 . (3.17)

(9)

3.4. Unconfined swimming

Once the shape of the swimmer is reconstructed as a function of time, the swimming velocity can be obtained by solving the Stokes equation outside of the swimmer. The velocity at the membrane is obtained as

u(α, t) = vs+dr(α, t)

dt , (3.18)

where v^s is the swimming velocity to be solved for. The deformation rate can be expressed using the complex notation as

dˆr(α, t)

dt =X

k6=0

drk

dt e^ik^α. (3.19)

The solution of the Stokes equation (2.2) can be expanded in multipoles as follows.

The complex coordinate is denoted as z = x + iy and its conjugate as z^∗. A general solution of (2.2) can be written (as recalled in appendix A) as

u = Aˆ (z) − z[A⁰(z)]^∗+ [B(z)]^∗, (3.20) where A(z) and B(z) are complex functions that are analytic in the fluid domain.

The functions A(z) and B(z) can be expanded in a series at infinity for the case of unconfined swimmer

A(z) =

∞

X

k=1

ak

z^k, B(z) =

∞

X

k=1

bk

z^k, (3.21a,b)

where b₁=0 because the volume of the swimmer is conserved (Re b₁=0) and the total torque acting on the swimmer is zero (Im b1=0). The expansion (3.21) incorporates the condition of the swimmer being force free because terms proportional to ln(z) would have to be added to functions A and B otherwise. The symmetry of the problem dictates that ak and bk be real. The coefficients ak and bk are expanded in powers of Γ as

ak=

∞

X

l=1

ak,lΓ^l^/2, b^k=

∞

X

l=1

bk,lΓ^l^/2. (3.22a,b) The next step is to substitute the expansions (3.21) and (3.22) into (3.20) and match the resulting flow field to the membrane velocity for z = ˆr(α).

Now we use the expression of the velocity field (3.20) evaluated at the swimmer surface and impose (3.18) (continuity of the velocity field). The shape coefficients rk

in (3.1) have been already determined (see (3.9)–(3.10)) as functions ofφk (defined in (3.6)). Expansion of (3.20) at the swimmer surface in powers of the small amplitude φk will generate various harmonics e^im^α (with m an integer). Using (3.18) and equating Fourier coefficients provides us relations between different coefficients. The zeroth Fourier harmonic relation yields the swimming velocity v^s. Note that the leading term of the swimming velocity requires expanding (3.20) to the order of Γ , and we obtain v^s= − π²cos(ωt)Γ¹^/2

L²₀[18 + 30 cos²(ωt)]¹^/2(a²+b2) + 3π³sin(ωt)Γ¹^/2

8L³₀[18 + 30 cos²(ωt)]¹^/2(a³+b3) + O(Γ³^/2).

(3.23)

(10)

Exploiting the relations following from the kth and −kth Fourier harmonics of (3.20) allows us to determine ak and bk at the desired order. The leading-order calculation yields

a2= 12L³₀ω cos(ωt)Γ¹^/2

π³[18 + 30 cos²(ωt)]³^/2 +O(Γ³^/2), (3.24) b₂= −L³₀ω[9 − 73 cos²(ωt)]Γ

8π³[3 + 5 cos²(ωt)]² +O(Γ³^/2), (3.25) b3= − 9L⁴₀ω sin(ωt)Γ¹^/2

π⁴[18 + 30 cos²(ωt)]³^/2 +O(Γ ), a3=O(Γ ). (3.26) Substituting (3.24)–(3.26) into (3.23) yields

vs=L0ωΓ [7 cos²(ωt) + 9]

[18 + 30 cos²(ωt)]² +O(Γ³^/2). (3.27) 3.5. Confinement effects

Equation (3.20) remains valid in the confined case but the domain where functions A(z) and B(z) can be considered analytic is now contained inside the strip interior to the walls. Therefore, their Laurent series about the origin contains all integer powers of z,

A(z) =

∞

X

k=1

ak

z^k +

∞

X

k=0

ckz^k, B(z) =

∞

X

k=1

bk

z^k +

∞

X

k=0

dkz^k. (3.28a,b) The coefficient c0 is assumed to be equal to zero below, which we can do without loss of generality: indeed, both c0 and d0 define a constant contribution to the velocity field when (3.28) is substituted into (3.20). The coefficients ck and dk can be expressed as a linear combination of ak⁰, a^∗_k0, bk⁰ and b^∗_k0 by imposing zero flow at the walls. These linear combinations come in form of infinite sums, the explicit analytical expression of which seems to be unavailable. We therefore resorted to approximate numerical calculations for this part. We simplify the following expressions by using the y → −y symmetry of the problem, which dictates the coefficients in expansion (3.28) to be real. In this case, the coefficient c₁ can be taken equal to zero: indeed, substituting expansions (3.28) into (3.20) shows that only the imaginary part of c1 enters the velocity field. The practical method can be described by the following algorithm:

(i) For each value of k, we calculate the residual velocity at the walls for a point singularity A(z) = 1/z^k, B(z) = 0 or A(z) = 0, B(z) = 1/z^k.

(ii) We then find the wall forces which would cancel the wall velocity computed in the previous step.

(iii) The final step is to calculate at the origin the velocity and its derivatives due to the wall forces found in the previous step.

The second and third steps have to be performed numerically, which can be done with high precision using Fourier representation of the residual velocity and forces at the wall. The derivatives of the velocity field define the coefficients ck and dk in (3.28). We can thus write using the linearity of (3.28)

ck=

∞

X

k⁰=1

c^ak_k⁰a_k⁰+c^bk_k⁰b_k⁰

W^k+k⁰ , dk=

∞

X

k⁰=1

d^ak_k⁰a_k⁰+d^bk_k⁰b_k⁰

W^k+k⁰ , (3.29)

(11)

where the coefficients c^ak_k⁰ and so on are dimensionless numbers obtained numerically as described above. These numbers do not depend on W because the width of the channel is the only length scale relevant in expansion (3.29). Expansion (3.29) corresponds to the no-slip condition on the channel walls. The form of the solution (3.29) can be understood simply. Indeed, when applying the no-slip boundary condition at the wall on (3.28), one obtains an infinite set of linear equations for the coefficients a_k, b_k, c_k and d_k, so that c_k and d_k are obviously linearly related to a_k, b_k. The dependence on W can be understood by rescaling z by W in (3.28), giving rise to W^−k and W^k in the first and second sums respectively.

The system is closed by imposing the velocity field (3.20) to be equal to the membrane velocity (3.18) at the boundary of the swimmer, as was done for the unconfined case. This allows us to compute the swimming velocity as a function of confinement.

As a practical demonstration of the method, we present here the intermediate steps for the derivation of the leading-order effect of the confinement. Retaining the terms up to O(C²n), we can write

d0=d^a2₀ a2+d^b2₀ b2

W² , d1=d^a1₁ a1

W² . (3.30a,b)

Here, we use that c1=0, c0=0 and b1=0, as discussed above. We also use that the coefficient d₀â1 (defining the contribution of a1 to d0) is equal to 0, as follows from the x → −x symmetry of the problem. Numerical calculation shows that d₀â2=3.7719, d^b2₀ = −6.8975 and dâ11 = −10.6695. Interestingly, we observe that

d₁^a1=d^b2₀ −d₀^a2. (3.31) This is not a numerical coincidence but a consequence of the invariance of the problem under translation along the x axis. Indeed, substituting z → z + dx, where dx is real, in (3.20) yields a velocity field that still satisfies the no-slip boundary conditions at the walls. Taking the derivative with respect to dx and setting dx = 0, yields that, if velocity field (3.20) satisfies no-slip boundary conditions at the channel walls for certain functions A(z) and B(z), then it would also satisfy them if we replace A(z) with A⁰(z) and B(z) with B⁰(z) − A⁰(z). This implies that expansions (3.29) remain valid if we replace ak+1 with −kak, ck−1 with kck, bk+1 with k(ak−bk) and dk−1 with k(dk−c_k). Applying this transformation to the first equation of (3.30) and comparing the result with the second one yields the relation (3.31).

The expression of the swimming velocity in the confined case reads v^s = − π²cos(ωt)Γ¹^/2

L²₀[18 + 30 cos²(ωt)]¹^/2(a²+b2) + 3π³sin(ωt)Γ¹^/2

8L³₀[18 + 30 cos²(ωt)]¹^/2(a³+b3) +d^a2₀ a2+d₀^b2b2

W² +O(C³n) + O(Γ³^/2). (3.32)

The coefficients ak and bk differ from the ones listed in (3.24)–(3.26) only in higher- order terms (i.e. to the desired order those expressions remain still valid here), except for a1 and b3, which read

a1= − 27ω sin(ωt)L²0Γ¹^/2 π²[18 + 30 cos²(ωt)]

1 − d₁^a1 L²₀

4π²W² +O(C³n)

+O(Γ ), (3.33) b3= − 9L⁴₀ω sin(ωt)Γ¹^/2

π⁴[18 + 30 cos²(ωt)]³^/2 − 27d^a1₁ L⁶₀ω sin(ωt)Γ¹^/2

16π⁶[18 + 30 cos²(ωt)]³^/2W² +O(Γ ). (3.34)

(12)

0.035

(a) (b)

0.030 0.025 0.020 0.015

√s/(øR) 0.010 0.005

0 π/2 π 3π/2 2π 0 π/2 π 3π/2 2π

0.015

0.010

0.005

0

-0.005

˝ = 0.054

˝ = 0.026

˝ = 0.010

Cn = 0.4 Cn = 0.2 Cn = 0.1 Cn = 0.04

FIGURE 2. (Colour online) The velocity of a swimmer vs as a function of time. The results of the analytical expansions (shown by continuous curves) are compared with the numerical results (shown by symbols). (a) Fixed confinement Cn=0.1. (b) Fixed excess perimeter Γ = 0.010.

Substituting the expressions for a_k and b_k into (3.32) yields the final expression v^s=L0ωΓ [7 cos²(ωt) + 9]

[18 + 30 cos²(ωt)]² +L³₀ωΓ [(73 cos²(ωt) − 9)d0^b2−9 sin²(ωt)d1^a1]

8π³W²[3 + 5 cos²(ωt)]² +O(Γ³^/2).

(3.35) 4. Results

4.1. Instantaneous velocity

In the previous section we have given the contribution to the swimming velocity to the leading order in confinement, which turns out to be of order C²_n (or inversely proportional to W²). However, the procedure outlined above can be applied in order to expand the swimming velocity to any order in confinement and excess perimeter.

In order to dispose of a wider range of applicability of the expansion, we have systematically performed it by retaining terms up to order O(Γ³) and O(C¹⁸n ). We have also validated the correctness and the convergence of the analytical results by full numerical simulations. The comparison of swimming velocity obtained analytically and by numerical simulations is shown in figure 2 as a function of time. As can be seen, the agreement is quite satisfactory. Nevertheless, increasing Cn and Γ beyond the shown values leads to a noticeable deviation of the analytical results from full numerical simulations (not shown in figure 2).

4.2. Average velocity

The instantaneous swimming velocity can be integrated in time to obtain the displacement over one stroke cycle 1X. The resulting series is sufficiently simple to be analysed in detail. In the absence of walls, the displacement 1X reads

1X = RΓ [2.1376 − 2.1621Γ + 1.4287Γ²+O(Γ³)] (4.1) as a function of the excess perimeter.

The dependence of 1X on confinement can be analysed by truncating the expansion in powers ofΓ at the leading order. The swimming speed is obtained as 1X/T, where T =2π/ω is the cycle period,

vs = RΓ

T 2.1376 − 0.39748C²n+18.933Cn⁴−52.316C⁶n+151.86Cn⁸−448.34C¹⁰n

+1374.5C¹²n −4230.7Cn¹⁴+12851.C¹⁶n −38256.C¹⁸n +O(C²⁰n ) + O(Γ²). (4.2)

(13)

0.1 0.2

Cn Cn

0.3 0.4 0.5

(a)0.20 (b)

0.15

Îx/R 0.10

0.05

0

˝ = 0.054

˝ = 0.026

˝ = 0.010

0 0.05 0.10 0.15 0.20

0.02160

0.02155

0.02150

0.02145

0.02140

FIGURE 3. (Colour online) Displacement of the swimmer after one stroke cycle as a function of the confinement for several values of the reduced volume. Comparison of analytical expansions (solid curves) with numerical simulations (symbols). (b) Shows the magnified view of the case Γ = 0.010 highlighting the initial decrease of the swimming velocity with increasing confinement.

Figure3shows the displacement 1X as a function of confinement for several values of the excess perimeter.

Comparison of numerical and analytical results in figure3 shows that the expansion (4.2) converges provided Cn is small enough. The ratio of the subsequent coefficients in (4.2) is close to −3. This suggests that the region of convergence of this expansion is about C²_n< 1/3. Higher-order terms in Γ in the expansion (4.2), however, have a smaller radius of convergence. This explains why the analytical curves in figure3start to deviate from the numerical results at lower values of Cn as Γ is increased.

5. Swimming with fixed shape

5.1. Propulsion by tangential flow

There are several swimmers for which the activity can be modelled by tangential flow.

A prototypical example is the paramecium (Zhang et al. 2015). The corresponding model is known as a squirmer, introduced by Lighthill (1952) and reconsidered later by Blake (1971), and has now become a very popular model (Lauga & Powers 2009). In a completely different field, that of mammalian cells, it was recently suggested (Aoun et al. 2019) that the swimming of immune cells (T-lymphocyte) is actuated by the retrograde flow of actin, which is transmitted to the outside fluid by transmembrane proteins, such as integrins. It is therefore an interesting question to see how confinement would affect this mode of locomotion.

The simplest model of a swimmer that uses tangential flow for locomotion in a shape-preserving manner is that of a squirmer. Here we consider a circular squirmer located in the centre of a channel. The flow at the boundary of the swimmer is written as

u(α) = v^s+v^t(α)t(α) ≡ v^s+(v⁰t ·t)t, (5.1) in the laboratory frame. Here v^t(α) is the local amplitude of the tangential velocity in the swimmer frame. We take vt =vt⁰ sinα ≡ v⁰t · t (with v⁰t pointing along the swimming x-direction), which is often considered as a model for ciliated organisms.

The solution process is similar to the one used for amoeboid swimmers but is simpler because the shape of the swimmer is precisely circular and the flow field is prescribed

(14)

0.2 0.4 Cn

0.6 0.8 1.0

1.0 0.8 0.6 0.4 0.2

√s/√s∞

0

FIGURE 4. (Colour online) Swimming velocity of a squirmer vs = |vsx| as a function of confinement, non-dimensionalized by the swimming velocity in unbounded fluid.

Continuous lines are calculated by (5.2), symbols are numerical simulations.

and constant in time. We thus need to replace (3.18), used for amoeboid swimmers, with (5.1) and recalculate vs using the same procedure as before. The calculations yield

vsx= −v⁰t

₁

2−0.86219C²n+0.76452Cn⁴−0.97679C⁶n+1.7022C⁸n−3.7459Cn¹⁰

+9.2095Cn¹²−23.687C¹⁴n +62.038C¹⁶n −163.97Cn¹⁸+435.74C²⁰n +O(C²²n ) . (5.2) The results obtained by (5.2) and by direct numerical simulations are shown in figure 4. Here again, we see a good agreement between the numerical and analytical results up to a confinement of approximately Cn=0.6.

5.2. Propulsion by normal flow

Following Zhu et al. (2013), we consider another case where the velocity at the surface is prescribed to be normal. This means that the swimmer pumps fluid through its membrane, ingesting it at the front region and expelling it at the rear, and this is the source of its motion. This swimmer will be referred to as pumper.

The velocity at the pumper surface is prescribed as

u(α) = ˜v^s+vⁿ(α)n(α), (5.3)

where we used the tilde for the swimming velocity for reasons that will become clear below; vⁿ(α) is a function of α, that we shall write (as for the tangential flow) as vⁿ(α) = vn⁰cos(α) ≡ v⁰n·n (with v⁰n pointing along the swimming x-direction). We can alternatively rewrite (5.3) as

u(α) = ˜v^s+(v⁰n·n)n = ˜v^s+v⁰n−(v⁰n·t)t = ˜V^s+(˜v⁰t ·t)t, (5.4) with ˜v⁰t= −v⁰n and ˜V_s= ˜vs+v⁰n. This means that the pumper problem is equivalent to the squirmer one (see (5.1)). This automatically yields the following identity between the two swimming speeds

v^sx vt⁰

= V˜sx

v˜t⁰

=v˜^sx+vn⁰

−vn⁰

. (5.5)

(15)

Whence, since vsx/vtx and ˜vsx/vnx are both negative, we can write v˜^s

vn⁰

+v^s vt⁰

=1. (5.6)

If the swimming speed is rescaled by the speed obtained for the unconfined geometry (denoted as vs^∞=vt/2 and ˜v^∞s =vn/2 for squirmer and pumper, respectively) the above identity reads (see (5.2))

v˜^s v˜s^∞

+ v^s vs^∞

=2. (5.7)

In three dimensions one has

v˜s

v˜^∞s

+2 vs

vs^∞

=3. (5.8)

The difference in numerical prefactors between (5.7) and (5.8) follows from the fact that the unbounded squirmer velocity in three dimensions is equal to 2/3 and is equal to 1/3 for the pumper (Zhu et al. 2013). The swimming speeds of the squirmer and pumpers have been calculated numerically by Zhu et al. (2013), and we have checked that identity (5.7) is consistent with their numerical results (see their figures 4 and 21).

We highlight the link between the squirmer and pumper by drawing the flow lines for each swimmer. The flow lines in a channel were calculated from full numerical solution and are shown in figure 5. The flow is exactly the same for both swimmer types in the wall frame of reference (figure 5a). The difference of the flow lines in the swimmer frame (figure 5b,c) arises from the difference of the swimming speed for squirmer and pumper.

6. Discussion

The main result of our analysis is (4.2) and (5.2) which presents the speeds of the amoeboid and squirmer swimmers as a function of the confinement. Interestingly, the only coefficient in (4.2) (amoeboid swimmer) that does not follow the same trend in term of amplitude is the one corresponding to C_n². This smallness is a consequence of the almost perfect cancellation of the two contributions which make up this coefficient:

the first one is due to the source dipole of the swimmer, which gets scattered off the walls, creating a flow opposite to the swimming direction. The second one is due to the stresslet of the swimmer, which creates a linear flow at the swimmer position when scattered off the walls. This linear flow modifies the swimming velocity by interacting with the deformations of the swimmer. The absence of this trend (smallness of the leading-order coefficient) for the squirmer is attributed to the absence of shape deformation.

The C²_n coefficient is negative for both swimmers, which means that vs decreases with increasing C_n for small values of C_n. However, since the coefficient of the C²n is anomalously small for an amoeboid swimmer, the next term becomes dominant for Cn much smaller than would otherwise be expected from the radius of convergence of (4.2). The C⁴_n coefficient in (4.2) is positive, which means that 1X increases with Cn for moderate confinements. We find that the minimum of 1X occurs at Cn = 0.1048. Numerical simulations (Wu et al. 2015, 2016) show that 1X goes through a maximum for high values of Cn and starts decreasing with increasing Cn. This region of parameters, however, is beyond the radius of convergence of (4.2)

(16)

4 3 2 1 0 -1

x/R y/R

-2 -3 -4 -1.5

1.00.5 0 -0.5 -1.0 y/R

-1.5 1.0 0.5 0 -0.5 -1.0 y/R (a) 1.5

(b)

(c) 1.0 0.5 0 -0.5 -1.0

-1.5

0.8 0.6 0.4 0.2 0 1.00 0.75 0.50 0.25 0 1.6 1.2 0.8 0.4 0

FIGURE 5. (Colour online) Flow lines around a confined circular swimmer. (a) Flow lines in the wall frame both for squirmer and pumper. (b) Flow lines in the swimmer frame for a squirmer. (c) Flow lines in the swimmer frame for a pumper. Numerical simulations.

Swimmers are located in the centre of the channel with W = 3R.

and thus cannot be analysed with the present calculation. Furthermore, the numerical simulations (Wu et al. 2015, 2016) were performed for large but finite values of the saturation number S. The approximation of fully saturated strokes becomes inadequate for strongly confined swimmers because of the slowdown of the lubrication dynamics as the gaps between the swimmer and the walls are narrow enough. Eventually, the channel becomes too narrow for the deformation amplitudes attained in fully saturated regime to be possible to deploy. The situation is different for the squirmer (or pumper) since the leading term in C²_n dominates causing the squirmer (pumper) to be monotonous.

The fact that the leading-order term is negative for the amoeboid swimmer and squirmer can be intuitively understood as follows. In order to move forward, the swimmer has to drag fluid from front to rear, and the presence of walls makes this operation more difficult, leading to a decrease of speed. The amoeboid swimmer has the ability to adapt its shape to reflection of the flow on the wall, taking advantage of its deformation to increase the gap between its shape and the wall. We believe that this is the source of the smallness of the C²_n coefficient. The pumper pumps fluid from the front to the rear so that it can drag the fluid through its body rather than on the side. The existence of a pumping alternative may be the reason for the increase of its velocity with confinement.

The present study uses a two-dimensional model of the swimmer. Nevertheless, some of the results can be expected to remain valid in three dimensions. The hydrodynamic interactions have a faster decay rate in three dimensions than in two dimensions: the flow field due to a point force depends as 1/r in three dimensions and as ln r in two dimensions on the distance r from the origin of the force. Thus, the flow field due to a force dipole decays as 1/r² in three dimensions, while the

(17)

flow field of a force quadrupole (equivalent to source dipole) decays as 1/r³. This suggests that the leading effect of confinement would be characterized by a term of order O(C³n) both for a tube and for a slit geometry.

Similarly to the two-dimensional case, the leading-order contribution of the confinement to the swimming velocity comes from two sources: the effect of the source dipole of the swimmer scattered off the walls and the linear flow at the swimmer position, which is created by the stresslet of the swimmer scattered off the walls. The coefficients with which these two effects enter the swimming velocity depend on the geometry of the confinement. A detailed analysis in three dimensions represents an important question for further research.

Acknowledgements

We thank CNES (Centre National d’Etudes Spatiales), ESA (European Space Agency), and the French–German university program ‘Living Fluids’ (grant CFDA- Q1-14). The simulations were performed on the Cactus cluster of the CIMENT infrastructure, which is supported by the Rhône-Alpes region (GRANT CPER07_13 CIRA).

Appendix A. Flow representation

Let us introduce a potential functionΛ(z, z^∗, t) to be related to the complex velocity field ˆu. The divergence-free condition ∇ · u = 0 can be satisfied by introducing a real- valued scalar potential Λ(z, z^∗, t), related to the velocity field ˆu as

u =ˆ ∂Λ

∂y −i∂Λ

∂x. (A 1)

Taking the curl of the first equation (2.2) entails that Λ is a biharmonic function

∇⁴Λ = 0. (A 2)

Upon using the relations between x, y and z, z^∗ one obtains that this equation transforms into

∂⁴Λ

∂z²∂z^∗2 =0. (A 3)

The general solution of which is

2Λ(z, z^∗) = f (z) + f1(z^∗) + z^∗g(z) + zg1(z^∗), (A 4) where the four functions on the right-hand side are arbitrary analytic functions of their respective complex argument (the factor 2 on the left-hand side is introduced for practical purposes). It is necessary that g1(z^∗) = (g(z))^∗ and f1(z^∗) = ( f (z))^∗ in order to guarantee that Λ is a real-valued function. This allows us to write finally that

Λ(z, z^∗) = Re[ f (z) + z^∗g(z)]. (A 5) The velocity field can then be written as

u = −2iˆ ∂Λ

∂z^∗ = −2i[g(z) + f (z)^∗0+zg(z)^∗0] ≡A(z) − zA⁰(z)^∗+B(z)^∗, (A 6) where A and B are two analytic functions obviously related to f and g.

REFERENCES

ACEMOGLU, A. & YESILYURT, S. 2014 Effects of geometric parameters on swimming of micro organisms with single helical flagellum in circular channels. Biophys. J. 106 (7), 1537–1547.