Stokesian swimming of a sphere at low Reynolds number

# Stokesian swimming of a sphere at low Reynolds number

## Abstract

Explicit expressions are derived for the matrices determining the mean translational and rotational swimming velocities and the mean rate of dissipation for Stokesian swimming at low Reynolds number of a distorting sphere in a viscous incompressible fluid. As an application an efficient helical propeller-type stroke is found and its properties are calculated.

###### pacs:
47.15.G-, 47.63.mf, 47.63.Gd, 87.17.Jj

## I Introduction

Stokesian swimming of a planar sheet by means of a transverse surface wave was first studied by Taylor (1). His analysis was soon extended to a squirming sphere by Lighthill (2). Blake corrected this work and applied it in a spherical envelope approach to ciliary propulsion (3). The calculations of mean translational swimming velocity and rate of dissipation are complicated, and for simplicity the analysis was restricted to axial strokes. Recently Pedley et al. (4) have included an axisymmetric azimuthal flow and showed that this can lead to a mean rotational swimming velocity. Shapere and Wilczek (5) derived a general expression for the mean translational swimming velocity of a sphere, but it is not easily applied, and analogous expressions for the mean rotational swimming velocity and the mean rate of dissipation are missing.

In earlier work (6) we have derived general expressions for mean translational swimming velocity and mean rate of dissipation for small amplitude swimming of a sphere in a fluid with inertia, but applied the analysis only to flow restricted to potential modes. Later we included viscous modes (7), and studied in particular the efficiency of swimming. Shapere and Wilczek (8) showed that a convenient definition of efficiency leads to the formulation of an eigenvalue problem involving two matrices, one for the mean translational swimming velocity and one for the mean rate of dissipation. We derived explicit expressions for the two matrices in a convenient representation, but for simplicity restricted the analysis to axial flows, in the manner of Lighthill (2) and Blake (3).

In the following we extend the analysis to the complete set of solutions of the steady state Stokes equations and derive explicit expressions for the integrals entering the matrices for mean translational and rotational swimming velocity and mean rate of dissipation. As an application we find an efficient helical propeller-type stroke and calculate its properties.

## Ii Flow equations

In this section we recall the basic equations of the bilinear theory of swimming. We consider a sphere of radius immersed in a viscous incompressible fluid of shear viscosity . At low Reynolds number and on a slow time scale the flow velocity and the pressure satisfy the Stokes equations

 η∇2\boldmathv−∇p=0,∇⋅\boldmathv=0. (1)

The fluid is set in motion by distortions of the spherical surface which are periodic in time and lead to swimming motion of the sphere as well as to a time-dependent flow field. The surface displacement is defined as the vector distance

 \boldmathξ=\boldmaths′−\boldmaths (2)

of a point on the displaced surface from the point on the sphere with surface . The fluid velocity is required to satisfy

 \boldmathv(\boldmaths+\boldmathξ(\boldmaths,t))=∂\boldmathξ(\boldmaths,t)∂t, (3)

corresponding to a no-slip boundary condition. The instantaneous translational swimming velocity , the rotational swimming velocity , and the flow pattern follow from the condition that no net force or torque is exerted on the fluid. We evaluate these quantities by a perturbation expansion in powers of the displacement .

To second order in the flow velocity and the swimming velocity take the form (9)

 \boldmathv(\boldmathr,t)=\boldmathv1(% \boldmathr,t)+\boldmathv2(\boldmathr,t)+...,\boldmathU(t)=\boldmathU2(t)+.... (4)

Both and vary harmonically with frequency , and can be expressed as

 \boldmathv1(\boldmathr,t) = \boldmathv1c(\boldmathr)cosωt+% \boldmathv1s(\boldmathr)sinωt, \boldmathξ(\boldmaths,t) = \boldmathξc(\boldmaths)cosωt+% \boldmathξs(\boldmaths)sinωt. (5)

Expanding the no-slip condition Eq. (2.3) to second order we find for the flow velocity at the surface

 \boldmathu1S(θ,φ,t) = \boldmathv1∣∣r=a=∂% \boldmathξ(θ,φ,t)∂t, \boldmathu2S(θ,φ,t) = \boldmathv2∣∣r=a=−\boldmathξ⋅∇\boldmathv1∣∣r=a, (6)

in spherical coordinates . In complex notation with the mean second order surface velocity is given by

 ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯\boldmathu2S(\boldmaths)=−12Re(\boldmathξ∗ω⋅∇)\boldmathvω∣∣r=a, (7)

where the overhead bar indicates a time-average over a period .

We consider periodic displacements such that the body swims in the direction. The time-averaged translational swimming velocity is given by (6)

 ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯\boldmathU2=¯¯¯¯¯¯U2\boldmathez,¯¯¯¯¯¯U2=−14π∫¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯\boldmathu2S⋅\boldmathezdΩ, (8)

where the integral is over spherical angles . Similarly the time-averaged rotational swimming velocity is given by (6)

 ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯\boldmathΩ2=¯¯¯¯¯¯Ω2\boldmath% ez,¯¯¯¯¯¯Ω2=−38πa∫(\boldmather×¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯\boldmathu2S)⋅\boldmathezdΩ. (9)

To second order the rate of dissipation is determined entirely by the first order solution. It may be expressed as a surface integral (6)

 D2=−∫r=a\boldmathv1⋅\boldmathσ% 1⋅\boldmatherdS, (10)

where is the first order stress tensor, given by

 \boldmathσ1=η(∇\boldmathv1+˜∇\boldmathv1)−p1\boldmathI. (11)

The rate of dissipation is positive and oscillates in time about a mean value. The mean rate of dissipation equals the power necessary to generate the motion.

## Iii Matrices

In our calculations it is convenient to expand the first order flow field and the pressure in terms of a basis set of complex solutions. The general solution of Eq. (2.1) which tends to zero at infinity and varies harmonically in time can be expressed as the complex flow velocity and pressure

 \boldmathvc1(\boldmathr,t) = −ωa∞∑l=1l∑m=−l[κlm\boldmathvlm(\boldmathr)+νlm\boldmathwlm(\boldmathr)+μlm\boldmathulm(\boldmathr)]e−iωt, pc1(\boldmathr,t) = −ωa∞∑l=1l∑m=−lκlmplm(\boldmathr)e−iωt, (12)

with complex coefficients and basic solutions (7),(10),(11)

 \boldmathvlm(\boldmathr) = (2l+2l(2l+1)^\boldmathAlm−2l−12l+1^\boldmathBlm)(ar)l, plm(\boldmathr) = 2η(2l−1)(−1)mPml(cosθ)eimφalrl+1, \boldmathwlm(\boldmathr) = −i^\boldmathClm(ar)l+1, \boldmathulm(\boldmathr) = −^\boldmathBlm(ar)l+2, (13)

with vector spherical harmonics in the notation of Ref. 11 (with in the normalization coefficient replaced by ), and with associated Legendre functions in the notation of Edmonds (12). The functions satisfy the Stokes equations (2.1), and the functions and satisfy these equations with vanishing pressure disturbance. For the solutions are axisymmetric and the functions are then identical with the functions introduced in Ref. 7. The solutions contain a factor , representing a running wave in the azimuthal direction for .

We consider a superposition of solutions of the form Eq. (3.1) with a single value of . The corresponding surface displacement takes the form

 \boldmathξ(\boldmaths,t)=−ia∞∑l=m[κlm\boldmathvlm(\boldmaths)+νlm% \boldmathwlm(\boldmaths)+μlm\boldmathulm(%\boldmath$s$)]e−iωt. (14)

Correspondingly we introduce the moment vector

 \boldmathψ=(κlm,νlm,μlm,κl+1,m,...). (15)

Then the mean second order swimming velocity is in the direction with value given by

 ¯¯¯¯¯¯U2=12ωa(\boldmathψ|B|%\boldmath$ψ$), (16)

with a dimensionless hermitian matrix . The mean second order rotational swimming velocity is in the direction with value given by

 ¯¯¯¯¯¯Ω2=34ω(\boldmathψ|C|\boldmathψ), (17)

with a dimensionless hermitian matrix . The time-averaged rate of dissipation can be expressed as

 ¯¯¯¯¯¯D2=8πηω2a3(\boldmathψ|A|\boldmathψ), (18)

with a dimensionless hermitian matrix . The matrix elements of the three matrices can be evaluated from Eqs. (2.6)-(2.11).

For we must put the moments equal to zero, because of the requirement that the swimmer exert no net force or torque on the fluid. Correspondingly for the matrices can be truncated by deleting the first two rows and columns. In the following we assume that this truncation has been performed.

The matrix turns out to be diagonal in the subscripts , as given by a factor . The basis elements will be indicated by a discrete index taking the values . Then the matrix at position has elements of the matrix of the form

 Alm=⎛⎜ ⎜ ⎜⎝alm000alm020alm110alm200alm22⎞⎟ ⎟ ⎟⎠. (19)

The nonvanishing elements are determined by the integrals

 ∫r=a\boldmathu∗km⋅(∇\boldmath% ulm)⋅\boldmatherdS = −4πafkm(k+1)(k+2)δkl, ∫r=a\boldmathu∗km⋅(∇\boldmath% vlm+˜∇\boldmathvlm−plm/η)⋅% \boldmatherdS = −8πafkm(k+1)(k+2)(2k−1)2k+1δkl, ∫r=a\boldmathv∗km⋅(∇\boldmath% ulm)⋅\boldmatherdS = −4πafkm(k+1)(k+2)(2k−1)2k+1δkl, ∫r=a\boldmathv∗km⋅(∇\boldmath% vlm+˜∇\boldmathvlm−plm/η)⋅% \boldmatherdS = −8πafkm(k+1)(2k3+k2−2k+2)k(2k+1)δkl,

with the factor

 flm=(l+m)!(l−m)!. (21)

For the expressions reduce to those given in Eq. (7.14) of Ref. 7. The integral is given by

 ∫r=a\boldmathw∗km⋅(∇\boldmathwlm+˜∇\boldmathwlm)⋅\boldmatherdS=−4πafkmk(k+1)(k+2)2k+1δkl. (22)

The and -matrices have matrices along the diagonal and are of the form

 Blm=⎛⎝B(1)lmB(2)lmB(3)lmB(4)lm⎞⎠, (23)

with matrices with the relation . The and -matrices have a checkerboard form with zeros on alternate positions. The integrals determining the -matrix are similar to those given in Eq. (7.7) of Ref. 7. Explicitly we find

 ∫r=a\boldmathv∗km⋅(∇\boldmath% vlm)⋅\boldmathezdS= − 4πa[gkm2(k+2)(2k−1)(2k+1)(2k+3)δk,l−1+glm(l+2)(2l−1)2(2l+1)(2l+3)δk,l+1], ∫r=a\boldmathv∗km⋅(∇\boldmath% wlm)⋅\boldmathezdS= 4πamfkm(k+1)(2k−1)2k+1δkl, ∫r=a\boldmathv∗km⋅(∇\boldmath% ulm)⋅\boldmathezdS= − 4πaglm(l+2)(2l+1)2l+3δk,l+1, ∫r=a\boldmathu∗km⋅(∇\boldmath% vlm)⋅\boldmathezdS= − 4πa[gkm2(k+2)2k+3δk,l−1+glm(l+2)(2l−1)2l+3δk,l+1], ∫r=a\boldmathu∗km⋅(∇\boldmath% wlm)⋅\boldmathezdS = 4πamfkm(k+1)δkl, ∫r=a\boldmathu∗km⋅(∇\boldmath% ulm)⋅\boldmathezdS= − 4πaglm(l+2)δk,l+1, (24)

with factor

 glm=(l+m+1)!(l−m)!. (25)

The integrals determining the -matrix are given by

 ∫r=a[\boldmather×(\boldmathv∗km⋅∇\boldmathvlm)]⋅\boldmathezdS = 4πiamfkm(k−2)(2k3+k2−2k+2)k2(2k+1)δkl, ∫r=a[\boldmather×(\boldmathv∗km⋅∇\boldmathwlm)]⋅\boldmathezdS = 4πia[gkm(k+2)2(2k−1)(2k+1)(2k+3)δk,l−1 − glml(l+2)(2l2+l+1)(l+1)(2l+1)(2l+3)δk,l+1], ∫r=a[\boldmather×(\boldmathv∗km⋅∇\boldmathulm)]⋅\boldmathezdS = 4πiamfkm(k+2)(2k−1)2k+1δkl, ∫r=a[\boldmather×(\boldmathw∗km⋅∇\boldmathvlm)]⋅\boldmathezdS = −4πia[gkm2k(k+2)(k+1)(2k+1)(2k+3)δk,l−1 − glm(l+2)(2l−1)(2l+1)(2l+3)δk,l+1], ∫r=a[\boldmather×(\boldmathw∗km⋅∇\boldmathwlm)]⋅\boldmathezdS = −4πiamfkm12k+1δkl, ∫r=a[\boldmather×(\boldmathw∗km⋅∇\boldmathulm)]⋅\boldmathezdS = 4πiaglml+22l+3δk,l+1, ∫r=a[\boldmather×(\boldmathu∗km⋅∇\boldmathvlm)]⋅\boldmathezdS = 4πiamfkm(k2−2)(2k−1)k(2k+1)δkl, ∫r=a[\boldmather×(\boldmathu∗km⋅∇\boldmathwlm)]⋅\boldmathezdS = 4πia[gkm(k+2)22k+3δk,l−1−glml(l+2)2l+3δk,l+1], ∫r=a[\boldmather×(\boldmathu∗km⋅∇\boldmathulm)]⋅\boldmathezdS = 4πiamfkm(k+2)δkl. (26)

The matrices and are real and symmetric. The matrix is pure imaginary and antisymmetric. In our earlier work (7) we have given explicit expressions for the upper left hand corners of the matrices and for , before deletion of the first two rows and columns. The explicit form of the truncated low order matrices for is given in an example of swimming by helical wave (13).

## Iv Eigenvalue problem

It is of interest to optimize the mean translational swimming velocity for given mean rate of dissipation . The optimization leads to a generalized eigenvalue problem of the form

 B|\boldmathψλ)=λA|\boldmath% ψλ). (27)

The maximum eigenvalue determines the optimum translational swimming velocity. The corresponding mean rotational swimming velocity can be found from the eigenvector by use of Eq. (3.6). In earlier work (7) we have shown that the maximum eigenvalue for is given by , but this obtains in the limit where moments of arbitrarily high number are involved. The absolute value of the elements of the corresponding eigenvector tends to a constant for large . Presumably the corresponding fine detail of the surface displacement is not physically relevant. In the following we consider finite-dimensional moment vectors, with at most equal to a bounded value . The corresponding matrices are -dimensional for . For the matrices are -dimensional. We denote the maximum eigenvalue with maximum equal to as .

In Fig. 1 we plot for as a function of . The maximum eigenvalue is maximum at , corresponding to axial symmetry. We find , not much less than . For the value is , and for it is .

In the eigenvector at the elements vanish, so that in optimum swimming the rotational modes with swirl are absent. The low order truncated matrix corresponding to reads explicitly

 C20=125⎛⎜ ⎜ ⎜⎝0010000010000000⎞⎟ ⎟ ⎟⎠. (28)

This shows that a moment vector with nonvanishing and elements leads to a nonvanishing . The elements correspond to the integrals for the pairs and in Eq. (3.15). Clearly there are similar couplings for . These couplings explain the nonvanishing in the study of Pedley et al. (4) for the axisymmetric case .

In earlier work (7) we have shown the absolute values of the moments of the optimal eigenvector for . We also showed the end of the displacement vector at various angles in the meridional plane . In Fig. 2 we show the real and imaginary parts of the radial displacement , and in Fig. 3 we show the real and imaginary parts of the component for the optimal eigenvector for . The eigenvector has been normalized to unity and its first component is chosen to be positive. The tangential component of is about twice as large as the radial component.

For the rotational swimming velocity does not vanish for the eigenvector corresponding to . In Fig. 4 we plot the ratio

 ρLm=(\boldmathξLm|CLm|\boldmathξLm)(\boldmathξLm|ALm|\boldmathξLm) (29)

as a function of for . This shows the mean rate of rotation for optimal strokes with maximum number equal to and different , normalized to equal power. In particular, at the ratio is , and at the ratio is . At the eigenvector is with nonvanishing moments .

The last example corresponds to a propeller-type first order flow. We consider more generally the case . The matrices are 3-dimensional and the eigenvalue equation can be solved in analytic form. The matrix is given by

 ALL=(L+1)(2L−1)!2L+1⎛⎜ ⎜⎝2L3+L2−2L+20L(L+2)(2L−1)012L2(L+2)0L(L+2)(2L−1)0L(L+2)(2L+1)⎞⎟ ⎟⎠. (30)

The matrix is given by

 BLL=iL(L+1)(2L)!2(2L+1)⎛⎜⎝0−2L+102L−102L+10−2L−10⎞⎟⎠. (31)

The matrix is given by

 CLL=(2L)!L(2L+1)⎛⎜⎝(L−2)(2L3+L2−2L+2)0L(2L3+L2−3L+1)0−L20L(2L3+L2−3L+1)0L2(L+2)(2L+1)⎞⎟⎠. (32)

We find for the maximum eigenvalue of the problem Eq. (4.1)

 λLL=√2L(2L+1)L+2. (33)

This tends to as , indicating an efficient swimmer. The corresponding eigenvector is

 \boldmathξLL=(0,1,−iL√2L(2L+1)), (34)

and the ratio is

 ρLL=L2+2L−2L2+3L+2. (35)

We denote the reduced power of the optimal swimmer with moments as ,

 PLL=(\boldmathξLL|ALL|\boldmathξLL)=L(L+1)(L+2)4L+2(2L)!. (36)

The time the swimmer needs to move over a distance equal to one diameter is (14)

 tLL=2a¯¯¯¯U2=4ωλLLPLL. (37)

During this time the swimmer rotates over the angle

 ¯¯¯¯Ω2tLL=3ρLLλLL. (38)

This is independent of the power, tends to as , and has a maximum at .

In Fig. 5 we show the imaginary part of the radial displacement for the eigenvector for . The real part vanishes. In Fig. 6 we show the real and imaginary parts of the component . In Fig. 7 we show the real and imaginary parts of the component . It is striking that these components have a simple dependence on the polar angle , concentrating the displacement along the equator. In the azimuthal direction there is a running wave given by the factor with .

## V Discussion

The calculation provides the complete set of matrix elements required in the discussion of translational and rotational swimming of a distorting sphere immersed in a viscous incompressible fluid, with no-slip boundary condition on the surface of the body. The calculation is in the spirit of Lighthill (2) and Blake 3, but allows general surface distortion. The last example of Sec. IV with a propeller type azimuthal surface wave demonstrates the power of the method.

One can imagine a wide range of possible surface waves. The optimization of translational swimming velocity for given power leads to an eigenvalue problem allowing selection of the best stroke. Structural considerations can lead to constraints imposed on the nature of the stroke. For example, one can require the surface displacement to be tangential (13) or radial (14).

The spherical geometry has the advantage that the complete set of solutions of the Stokes equations is known in terms of elementary functions. The calculation can serve as a guide to the study of more complicated geometry, for example with the goal of finding the optimal swimming modes of a spheroid.

## Figure captions

### Fig. 1

Plot of the maximum eigenvalue as a function of . This characterizes the mean translational swimming velocity of the optimal swimmer for at each value of for given power.

### Fig. 2

Plot of the real (solid curve) and imaginary (dashed curve) parts of the radial displacement in the meridional plane at time as a function of polar angle for the optimal swimmer with .

### Fig. 3

Plot of the real (solid curve) and imaginary (dashed curve) parts of the tangential displacement in the meridional plane at time as a function of polar angle for the optimal swimmer with .

### Fig. 4

Plot of the ratio , defined in Eq. (4.3), for as a function of . This characterizes the rate of steady rotation of the optimal swimmer for at each value of for given power.

### Fig. 5

Plot of the radial displacement in the meridional plane at time as a function of polar angle for the optimal swimmer with .

### Fig. 6

Plot of the real (solid curve) and imaginary (dashed curve) parts of the tangential displacement in the meridional plane at time as a function of polar angle for the optimal swimmer with .

### Fig. 7

Plot of the real (solid curve) and imaginary (dashed curve) parts of the azimuthal displacement in the meridional plane at time as a function of polar angle for the optimal swimmer with .

### References

1. G. I. Taylor, Analysis of the swimming of microscopic organisms, Proc. R. Soc. London A 209, 447 (1951).
2. M. J. Lighthill, On the squirming motion of nearly spherical deformable bodies through liquids at very small Reynolds numbers, Comm. Pure Appl. Math. 5, 109 (1952).
3. J. R. Blake, A spherical envelope approach to ciliary propulsion, J. Fluid Mech. 49, 209 (1971).
4. T. J. Pedley, D. R. Brumley, and R. E. Goldstein, Squirmers with swirl - a model for Volvox swimming, arXiv:1512.02931[cond-mat.soft].
5. A. Shapere and F. Wilczek, Geometry of self-propulsion at low Reynolds number, J. Fluid Mech. 198, 557 (1989).
6. B. U. Felderhof and R. B. Jones, Small-amplitude swimming of a sphere, Physica A 202, 119 (1994).
7. B. U. Felderhof and R. B. Jones, Optimal translational swimming of a sphere at low Reynolds number, Phys. Rev. E 90, 023008 (2014).
8. A. Shapere and F. Wilczek, Efficiencies of self-propulsion at low Reynolds number, J. Fluid Mech. 198, 587 (1989).
9. B. U. Felderhof and R. B. Jones, Inertial effects in small-amplitude swimming of a finite body, Physica A 202, 94 (1994).
10. R. Schmitz and B. U. Felderhof, Creeping flow about a spherical particle, Physica A 113, 90 (1982).
11. B. Cichocki, B. U. Felderhof, and R. Schmitz, Hydrodynamic interactions between two spherical particles, PhysicoChem. Hyd. 10, 383 (1988).
12. A. R. Edmonds, Angular Momentum in Quantum Mechanics (Princeton University Press, Princeton, NJ, 1974).
13. B. U. Felderhof, Spinning swimming of Volvox by tangential helical wave, arXiv:1601.00755[physics.flu-dyn].
14. B. U. Felderhof, Stokesian swimming of a sphere by radial helical surface wave, arXiv:1601.03151[physics.flu-dyn].
105030