PMCCPMCCPMCC

Search tips
Search criteria 

Advanced

 
Logo of procaThe Royal Society PublishingProceedings AAboutBrowse by SubjectAlertsFree Trial
 
Proc Math Phys Eng Sci. 2016 August; 472(2192): 20160425.
PMCID: PMC5014114

Optimization and small-amplitude analysis of Purcell's three-link microswimmer model

Abstract

This work studies the motion of Purcell's three-link microswimmer in viscous flow, by using perturbation expansion of its dynamics under small-amplitude strokes. Explicit leading-order expressions and next-order correction terms for the displacement of the swimmer are obtained for the cases of a square or circular gait in the plane of joint angles. The correction terms demonstrate the reversal in movement direction for large stroke amplitudes, which has previously only been shown numerically. In addition, asymptotic expressions for Lighthill's energetic efficiency are obtained for both gaits. These approximations enable calculating optimal stroke amplitudes and swimmer's geometry (i.e. ratio of links’ lengths) for maximizing either net displacement or Lighthill's efficiency.

Keywords: microswimmers dynamics, Purcell's three-link swimmer, perturbation analysis, optimization

1. Introduction

The study of micrometre-size swimmers dynamics has in recent years become a highly active research area and has possible implications on understanding the motion of swimming microorganisms and biological infections [1,2]. The considerable advances made in the field of micro- and nanotechnology have promoted the possibility of manufacturing miniature robotic devices that operate in these small scales [35]. Such mini-robots may have many applications in medicine, performing medical procedures in a minimally invasive way and delivering drugs with high precision [68]. All this requires an understanding of swimming dynamics at the low Reynolds number regime.

The Reynolds number represents the ratio of the inertial forces to the viscous ones. When dealing in microfluidics, where the characteristic lengths are extremely small, hydrodynamic forces are typically governed by very low Reynolds numbers (Re[double less-than sign]1) [9]. The result of this is that the strategy of motion in this regime needs to be drastically different from the familiar motion of larger organisms, such as fish, that rely on imparting momentum to the surrounding fluid. In his famous lecture [10], Purcell introduced the ‘Scallop theorem’, which states that a swimmer that changes its shape and then changes back to the original shape by reversing the same sequence will return to the point it started its motion with no regard to the speed in which any part of the motion is made. Any reciprocal motion of the swimmer will result in zero net translation. There are several ways of overcoming the scallop theorem and generating net motion in a low Reynolds regime. One of them is continuously performing a unidirectional rigid body motion such as a rotating corkscrew. In this way, there is no reciprocal motion and the result is generation of net motion of the swimmer. This is precisely the method used by the Escherichia coli bacteria to propel itself [11]. Another way for overcoming the scallop theorem is by making a non-reciprocal periodic shape change, which will be henceforth called a ‘gait’. The simplest version of such a swimmer, known as ‘Purcell's three-link swimmer’ (figure 1), was suggested by Purcell in [10] and can be seen as a simplified version of the travelling wave ‘Taylor sheet’ [12] in two dimensions which is discretized to have only 2 degrees of freedom. Purcell's swimmer comprises three rigid links connected by two rotary joints (figure 1). Purcell indicated that this swimmer could propel itself along a straight line by alternately rotating its front and back links in a non-reciprocal way (square gait, figure 2a). Through symmetry considerations alone, it can be shown that the three-link swimmer will move along the x-axis when using the shape changes suggested by Purcell. Purcell claimed that determining the direction of net motion (i.e. forward or backward) is trivial and left it as ‘an exercise for the student’. Only 26 years later, Becker et al. [13] obtained an explicit formulation for the dynamics of this microswimmer, and surprisingly found that the direction of net motion actually depends on the stroke amplitude of joint angles. Specifically, for small amplitudes the swimmer will move in one direction but for larger amplitudes the swimmer will move in the opposite direction. Additionally, Becker et al. [13] studied Lighthill's energetic efficiency, which is roughly equivalent to the maximal mean speed achievable under a given mean expenditure of mechanical power. It is worth noting that both [10,13] considered only the case where the joint angles rotate alternately, creating a square gait and did not study other possible periodic shape change, such as the circular gait (figure 2b), at which the two joint angles oscillate sinusoidally with a quarter-period phase shift.

Figure 1.
Purcell's three-link swimmer. (Online version in colour.)
Figure 2.
Gait representation in joint angle plane. (a) Square gait and (b) circular gait. (Online version in colour.)

Some previous works have examined Purcell's swimmer from different perspectives. In [14], Gutman and Or analyse the symmetries of Purcell's swimmer, derive conditions on gaits that result in movement along the swimmer's principal directions and present motion experiments of a macro-scale robotic swimmer in a highly viscous fluid. Avron & Raz [15] and Hatton and colleagues [16,17] use tools of differential geometry in order to obtain geometric visualization of performance measures in shape space (plane). By observing curvature maps of the swimmer's dynamics, those works can qualitatively explain the reversal in direction of motion for large-amplitude strokes, and also provide guidelines that help in obtaining optimal gaits. Tam & Hosoi [18], through numerical computation only, found gaits that achieve optimal translation during a period as well as energetically optimal gaits. They also found the optimal geometry of the swimmer, i.e. ratio between the swimmer's links lengths, for these two optimality criteria, again only through numeric calculations.

Aside from numerical computations, another useful approach is obtaining closed-form expressions under some simplifying scaling assumptions, by using asymptotic analysis. In microswimmers, this approach dates back to the analysis of Taylor [12], who showed that the swimming speed of Taylor's sheet scales quadratically with the wave's amplitude at leading order. In [19], this concept was employed on analysis of self-propulsion of spherical swimmers performing a squirming motion. For the shape-changing three-link swimmer, one has to use the method of perturbation expansion [20] under the assumption of small stroke amplitudes. This method was used in [13] for obtaining leading-order expressions for the motion of Purcell's swimmer and finding its optimal geometry. In [21], the leading-order expressions were derived using a similar approach and then used to find optimal geometry of the swimmer for the square gait and for derivation of an optimal polygonal gait under small joint angles. Interestingly, there are differences in the results found in [13,18,21] in terms of the swimmer's optimal geometry. Recently, perturbation expansion has become a more common tool for analysis of simple microswimmers motion in various cases such as a three-link swimmer with a passive elastic joint [22], a two-link magnetically actuated microswimmer [23,24] and wobbling of a magnetically actuated helix [25].

The goal of this work is to introduce a systematic method for analysing the dynamics of Purcell's three-link swimmer using perturbation expansion and to exploit this method to perform optimizations on the swimmer geometry and stroke amplitude. The main new contributions of this work compared with previous literature are (i) formulation of leading-order expression for the displacement under the circular gait (figure 2b); (ii) derivation of the next-order correction terms for the displacement under both square and circular gaits, which explicitly show the reversal in movement direction first found in [13], as well as obtaining an approximation for the optimal stroke amplitude; (iii) obtaining leading-order expressions for Lighthill's energetic efficiency of the two presented gaits and also a next-order term for the square gait, as well as using these expressions to perform optimization of the swimmer's geometry and stroke amplitude for maximal efficiency. (iv) Resolving the previously unexplained differences between the findings in [13,18,21] regarding the optimal geometry.

The paper is organized as follows. The next section presents Purcell's three-link swimmer model and its equations of motion, the two gaits in question and the concept of Lighthill's energetic efficiency. Section 3 focuses on perturbation expansion and presents a systematic way for deriving the net displacement of the swimmer in the form of a power series for any gait, and the first two terms are found explicitly for the square and circular ones. The optimal geometry and stroke amplitude for maximal displacement are also found. This section also includes the use of perturbation expansion to find the period time under constant joint torque for the square gait, and under constant power for both square and circular gaits. In §4, the previous results are used in order to obtain the leading-order terms of the energetic efficiency of the two gaits, as well as next-order term for the square gait. Optimal geometry and amplitude are then found for maximal efficiency. Finally, §5 offers our conclusions, indicating possible future extensions and consequences.

2. Problem formulation

In this section, we introduce Purcell's swimmer model, formulate its dynamic equations of motion, define the square and circular gaits and, finally, review the concept of Lighthill's energetic efficiency. The swimmer in question consists of three thin rigid links with lengths l0,l1,l2, with l1=l2, l=l0+l1+l2 and η=l0/l. The links are connected by two rotary joints whose angles are denoted by ϕ1 and ϕ2 (figure 1). The shape of the swimmer will be described by these two angles ϕ=(ϕ1,ϕ2)T. It is assumed that the swimmer's motion is confined to the xy plane. The planar position and orientation of the centre link are denoted by q=(x,y,θ)T. The swimmer is submerged in an unbounded fluid domain whose motion is governed by Stokes equations under low Reynolds number [9]. The velocity of the ith link is described by the linear velocity of its centre vi and the link's angular velocity ωi, which are augmented in the vector Vi=(vi,ωi)R3. Similarly, let Vb=(vb,ωb) denote the linear and angular velocities of a body-fixed reference frame, attached to the central link. In addition, we define the indicator matrix Iij, such that Iij=±1 if the location of the ith link with respect to the body frame is affected by the jth joint, with the sign determined by the direction of the joint axis. Iij=0, if the jth joint does not affect the ith link. The internal torques (moments) acting on the joints are denoted by τ1 and τ2.

(a) Formulation of equations of motion

The kinematic relation between body velocity, joints velocities and links velocities is given by

vi=vb+ωbz×(rir0)+jIijϕ˙jz×(ribj)andωi=ωb+jIijϕ˙j,}
2.1

where ri is the position of the centre of the ith link and bj is the position of the jth joint. In matrix form, the velocity Vi is related to the body velocity Vb and shape velocity ϕ˙ through

Vi=Ti(q,ϕ)q˙+Ei(q,ϕ)ϕ˙.
2.2

The matrices Ti(q,ϕ) and Ei(q,ϕ) for i=0,1,2 are given by

T0=[100010001],E0=[000000]T1=[100.5l0sinα00.5l1sinα1010.5l0cosα0+0.5l1cosα1001],E1=[0.5l1sinα100.5l1cosα1010]andT2=[100.5l0sinα0+0.5l2sinα2010.5l0cosα00.5l2cosα2001],E2=[00.5l2sinα200.5l2cosα201],}
2.3

where α0=θ, α1=ϕ1θ and α2=θ+ϕ2 are the absolute orientation angles of each link. Next, we invoke resistive force theory [26,27], which states that the viscous drag force fi and torque mi on the ith slender link under planar motion are proportional to its linear and angular velocities. Thus, we can write the expression for the drag force and torque exerted on each link:

fi=ctli(viti)ticnli(vini)niandmi=112cnli3ωi,}
2.4

where ti=(cosαi,sinαi)T is a unit vector in the axial direction of the ith link and ni=(sinαi,cosαi)T is a unit vector in the normal direction. The resistance coefficients for the normal and axial directions are cn=2ct=4πμ/log(lc/a), where μ is the dynamic viscosity of the fluid, a is the radius of the links and lc is a characteristic length. As the ratio of the link's length to its radius is assumed to be very large, the difference in the resistance coefficients between the links is very small and will be neglected. It is also assumed that effects of hydrodynamic interaction between the links are negligible. Denoting the vector of forces and torques on the ith link as Fi=(fi,mi), the relation (2.4) can be written in matrix form

Fi=Ri(q,ϕ)Vi,
2.5

where

Ri(q,ϕ)=ct(i)li[1+sin2αicosαisinαi0cosαisinαi1+cos2αi00016li2]
2.6

is called the resistance tensor. The total drag forces and torques acting on the swimmer's body are given by

fb=i=02fi,mb=i=02(mi+((rbri)×fi)z).
2.7

Using the matrices Ti from the kinematic relation (2.2) and augmenting in a vector Fb=(fb,mb), (2.7) is written in matrix form as

Fb=i=02TiTFi=i=02TiTRi(TiVb+Eiϕ˙).
2.8

Denoting

Rbb=i=02TiTRiTiandRbu=i=02TiTRiEi,
2.9

equation (2.8) becomes

Fb=RbbVb+Rbuϕ˙.
2.10

Our choice of coordinates of body location q implies that the body velocity satisfies q˙=Vb. Assuming quasi-static motion, the swimmer is in static equilibrium Fb=0. Substituting this into (2.10) then yields the nonlinear differential equations that govern the swimmer's dynamics:

q˙=G(q,ϕ)ϕ˙,
2.11

where

G(q,ϕ)=Rbb1Rbu.
2.12

As there are no external boundary conditions in the unbounded fluid domain, the swimmer's velocity expressed in body-fixed reference frame is independent of the position q and thus (2.11) can be written as (cf. [14])

q˙=D(θ)G(ϕ)ϕ˙,whereD(θ)=[cosθsinθ0sinθcosθ0001].}
2.13

Note that (2.13) implies that the angular velocity of the swimmer θ˙ is independent of θ. The matrix G(ϕ) obeys some symmetry relations due to the swimmer's structure (see [14] for details):

G(ϕ)=[100010001]G(ϕ))andG(Sϕ)=[100010001]G(ϕ)S,}
2.14

where S=[0110] represents interchanging between the joint angles ϕ1 and ϕ2.

Next, we compute the actuation torques acting at the joints, which are denoted by the vector τ=(τ1,τ2)T. Owing to static equilibrium of the partial kinematic chain ending at the jth joint, these torques are balanced by hydrodynamic forces and torques fi and mi, giving rise to the following relation:

τj=iIij(mi+((ribj)×fi)z).
2.15

This relation can be written in matrix form using the previously defined matrices Ei:

τ=i=02EiTFi.
2.16

Substituting (2.2), (2.5) and (2.11) in (2.16) gives

τ=i=02EiTRiVi=i=02EiTRi(Tiq˙i+Eiϕ˙i)=(i=02EiTRiEii=02EiTRiTi(i=02TiTRiTi)1i=02TiTRiEi)ϕ˙
2.17

denoting Ruu=i=02EiTRiEi; this equation can be written as

τ=(RuuRbuTRbb1Rbu)ϕ˙
2.18

with Rbb and Rbu as defined earlier in (2.9). Denoting

W(ϕ)=RuuRbuTRbb1Rbu,
2.19

we have

τ=W(ϕ)ϕ˙.
2.20

The resulting dynamic equations in (2.11) and (2.20) are strongly nonlinear. Nevertheless, they can be expanded by assuming small-amplitude changes of the angles about ϕ=0, as explained in §3.

(b) Periodic gaits

This work considers time-periodic inputs of shape changes ϕ(t), called gaits, which represent closed loops in the plane of joint angles. In particular, we focus on two specific possible gaits, square and circular (figure 2), that are presented here with ε as a scaling factor of the stroke amplitude, which will later be assumed small. The time function of the relative angles between the links can be written as ϕi(t)=εsi(t), where si represents the ‘unscaled’ shape trajectory. Let us also define the vector s=(s1,s2)T. For the cases of square and circular gaits, the joint angles can be written in the unscaled form as

square:s1(t)={1,t[0,2](3t),t[2,4]1,t[4,6](t7),t[6,8],s2(t)={(1t),t[0,2]1,t[2,4](t5),t[4,6]1,t[6,8]
2.21

and

circular:s1(t)=sin(t+π4),s2(t)=cos(t+π4),t[0,2π].
2.22

These equations describe circular and square-shaped trajectories with stroke amplitude of 1, which are then scaled to stroke of ε by setting ϕi(t)=εsi(t). As the equation of motion (2.11) is time invariant, the net motion is independent of time parametrization of the gait. Here, time parametrization of the gaits was chosen arbitrarily so that the period times are T=8 for the square gait and T=2π for the circular one. The net displacement in the x-direction of one full period will be denoted X and is calculated through X=0Tx˙dt. The mean swimming speed is denoted by V =X/T.

(c) Mechanical power and efficiency

The energetic efficiency of a stroke can be defined in a several different ways, cf. [28]. First, we write an expression for the power exerted by the swimmer. The mechanical power dissipated by the fluid's viscous drag forces and torques on all three links is

P=i=02FiTVi=i=02ViTRiVi.
2.23

On the other hand, the mechanical power expended internally by the actuation torques is

P=τTϕ˙=ϕ˙TWϕ˙.
2.24

These last two expressions are equivalent, which can be proved using the relations (2.19) and (2.20). The total work can be calculated from here by W=0TP(t)dt.

Owing to time-invariance of the swimmer's dynamics, a well-known observation (cf. [13,18]) states that the energy expenditure under a given gait trajectory can be made arbitrarily small by pacing along the trajectory in a sufficiently slow rate. Thus, energy per unit distance cannot serve as a reasonable performance measure. Following Becker et al. [13] and Tam & Hosoi [18], we use an energetic efficiency criterion similar to that defined by Lighthill in [29] which is the ratio of average power P¯=W/T exerted by the swimmer to the power needed to drag the swimmer as a rigid body at the same mean speed:

ξ~=ctlV2P¯=ctlX2P¯T2.
2.25

This definition is non-unique for a given gait trajectory ϕ(σ(t)) and by varying the gait's time parametrization σ(t) the average power also changes. Nevertheless, a known result from Becker et al. [13] proves that minimal average power P¯ for a given gait is obtained by choosing a time parametrization σ(t) for which the instantaneous power P(t) is kept constant. Using this particular time parametrization (which is unique up to multiplying by a positive factor) maximizes the energetic efficiency ξ~ for a given gait, and this efficiency is uniquely determined for any gait trajectory.

The period time T under a constant mechanical power P=Po can be found using the following derivation. Consider the gait ϕ(σ)=ϕ(σ(t)), where σ[set membership][σ0,σ1] is a geometric parameter along the gait's trajectory and σ(t) is its time parametrization. Using (2.24), the instantaneous power P(t) can be written as

P(σ(t))=(ϕ(σ)TW(σ)ϕ(σ))σ˙2,whereϕ(σ)=ϕσ
2.26

denoting F(σ)=ϕ′(σ)W(σ)ϕ′(σ) we have

σ˙=dσdt=PF(σ)dt=F(σ)Pdσ.
2.27

This transformation is well defined for any non-degenerate trajectory and time parametrizations such that ϕ′(σ)≠0 and dσ/dt>0, for all σ and t. The period time of a gait under constant mechanical power P=Po is thus obtained as

Tpo=1Poσ0σ1F(σ)dσ.
2.28

Substituting (2.28) into the expression for Lighthill's efficiency (2.25) under a constant power P¯=Po, it is clear that the Po cancels out and so the calculations can be done for Po=1 without loss of generality. Thus, Lighthill's efficiency of a gait ϕ(σ) is uniquely given as

ξ=X2Tp2,
2.29

where Tp now denotes the period time under constant power of Po=1, as given in (2.28), and the constants ctl from (2.25) are dropped.

3. Perturbation expansion of the dynamics

In order to find the leading-order expressions for the swimmer's motion using perturbation expansion [20], the dynamics of the swimmer are expanded as power series of the stroke amplitude ε. First, expansion of the swimmer's orientation angle θ(t) is obtained, followed by expansions of the instantaneous body velocity and of the resulting net displacement X under each specific gait.

(a) Expansion of swimmer displacement

The position and orientation of the swimmer are now expanded into a power series in ε as

q(t)=εq(1)(t)+ε2q(2)(t)+ε3q(3)(t)+.
3.1

And for each coordinate

x(t)=εx(1)(t)+ε2x(2)(t)+ε3x(3)(t)+,
3.2

y(t)=εy(1)(t)+ε2y(2)(t)+ε3y(3)(t)+
3.3

andθ(t)=εθ(1)(t)+ε2θ(2)(t)+ε3θ(3)(t)+.
3.4

From equation (2.13), it can be observed that θ˙ is independent of θ, and thus the ODE for θ can be written as

θ˙=j=12G3j(ϕ)ϕ˙j=εj=12G3j(ϕ)s˙j.
3.5

Using Taylor expansion we get

θ˙=j=12[G3j(0)+(ϕ1ϕ1|(0,0)+ϕ2ϕ2|(0,0))G3j+12!(ϕ1ϕ1|(0,0)+ϕ2ϕ2|(0,0))2G3j+]ϕ˙j.
3.6

Substituting the expansion for θ and ϕi(t)=εsi(t):

εθ˙(1)+ε2θ˙(2)+ε3θ˙(3)+=j=12[G3j(0)+ε(s1ϕ1|(0,0)+s2ϕ2|(0,0))G3j+ε212!(s1ϕ1|(0,0)+s2ϕ2|(0,0))2G3j+]εs˙j.
3.7

Now we can collect terms of different orders in ε.

θ˙(1)(t)=j=12G3j(0)s˙j(t),
3.8

θ˙(2)(t)=j=12(G3jϕ1(0)s1+G3jϕ2(0)s2)s˙j(t)
3.9

andθ˙(3)(t)=j=12(122G3jϕ12(0)s12+122G3jϕ22(0)s22+2G3jϕ1ϕ2(0)s1s2)s˙j(t).
3.10

Explicit expressions for the derivatives are given in the electronic supplementary material. The first-order derivatives at the origin are zero as G31, G32 are even functions in (ϕ1,ϕ2) and so θ˙(2)=0 in (3.9), implying that θ(2)(t)=0 is zero (under zero initial conditions). Now, using the expression for θ(t) we can calculate expansions of the net motion q(t). First, we expand the matrices in (2.13) as

D(θ)=I+θJ+12!θ2J2+13!θ3J3+andJ=[010100000]}
3.11

and

G(ϕ)=G(0)+ε(s1ϕ1|(0,0)+s2ϕ2|(0,0))G+ε212!(s1ϕ1|(0,0)+s2ϕ2|(0,0))2G+.
3.12

This, of course, is the same process done for the last row of G in (3.6). Expanding equation (2.13) by substituting expansions for θ(t), D(θ), G(ϕ) and also of q and ϕi(t)=εsi(t) and collecting orders of ε, one obtains

ε1:q˙(1)=G(0)s˙(t)ε2:q˙(2)=[θ(1)JG(0)+(s1ϕ1|(0,0)+s2ϕ2|(0,0))G]s˙(t)ε3:q˙(3)=[12θ(1)2J2G(0)+θ(1)J(s1ϕ1|(0,0)+s2ϕ2|(0,0))G+12(s1ϕ1|(0,0)+s2ϕ2|(0,0))2G]s˙(t)ε4:q˙(4)=[(θ(3)J+13!θ(1)3J3)G(0)+12θ(1)2J2(s1ϕ1|(0,0)+s2ϕ2|(0,0))G+12θ(1)J(s1ϕ1|(0,0)+s2ϕ2|(0,0))2+16(s1ϕ1|(0,0)+s2ϕ2|(0,0))3G]s˙(t).}
3.13

Explicit expressions for the derivatives of elements of G are given in the electronic supplementary material. Owing to symmetry of the swimmer and of the gaits, it is known that there will only be net translation along x-direction, which is the axis of the centre link, while rotation and translation in y-direction are cancelled out [14]. Thus, expansions of net displacement X for specific gaits are calculated next.

(b) Gait-specific expressions of the displacement

Once the series expansion for the instantaneous velocities of the swimmer is derived, it is possible to obtain the net displacement X over a period for a specific gait via integration. The following proposition summarizes the results for the square and circular gaits.

Proposition 3.1 —

For a symmetric, three linked ‘Purcell swimmer’ performing a square or circular gait with amplitude ε, the leading-order term and next-order correction for the displacement X over one full stroke in the direction of x axis are given as

X=f2(η)ε2f4(η)ε4+O(ε6).
3.14

With

f2(η)=C2ηl(1η)3(η+3)andf4(η)=C4ηl(1η)3(η7+3η610η522η4+29η3+95η2+44η+20),}
3.15

where for the square gait given in (2.21), we have C2=14,C4=1192 and for the circular gait given in (2.22) C2=π/16, C4=π/1024.

The proof of this proposition is given in the electronic supplementary material. It can easily be seen that the leading-order term of X is quadratic in ε and that the next-order corrections are of opposite sign to the leading-order terms. This implies that the displacement grows monotonically with ε for small stroke amplitudes, whereas for larger amplitudes with ε>1 the O(ε4)-term causes reversal in the direction of net motion. This, in turn, indicates the existence of a locally optimal value of ε that achieves maximal displacement. Note that, a leading-order term for displacement under the square gait is found in [21] using Lie brackets, and the results are identical to those given here.

As a demonstration of the utility of proposition 3.1, figure 3 shows plots of the displacement X versus stroke amplitude ε with link ratio of η=13 for both gaits (figure 3a, square, and and33b, circular). The solid curves are obtained from numerical integration of the nonlinear equations of motion (2.11), whereas the dashed and dash-dotted curves are O(ε2) and O(ε4) approximations, respectively. While the O(ε2)-approximation is monotonic in ε and works only for small stroke amplitudes, the O(ε4)-approximation with next-order correction captures the reversal in the direction of the displacement for intermediate amplitudes, though there is an increasing deviation from the numerical values for larger amplitudes, which is of O(ξ6).

Figure 3.
Displacement over one cycle X for η=13 as a function of stroke amplitude ε for (a) square gait and (b) circular gait. Solid curves, numerical integration. Dashed curves, O(ε2) approximation. Dash-dotted curves, O(ε2) approximation. ...

(c) Optimal geometry and stroke amplitude for maximal displacement

First, we discuss the dependence of net displacement X on the swimmer's geometric ratio η and derive its locally optimal values for both gaits, based on equations (3.14) and (3.15). As a demonstration, plots of X versus η under the square gait are shown in figure 4a,b for stroke amplitudes of ε=π/4 and ε=π/2, respectively. It can be seen that the O(ε2) approximation (dashed) has a larger deviation from the exact value obtained by numerical integration (solid curves), compared with that of the O(ε4) approximation (dash-dotted), where the deviations are further increased for larger amplitudes ε. Nevertheless, in all cases it is obvious that the displacement X vanishes at extreme cases of η→0 or η→1, where either the middle link or side links vanish, and that an optimal value of η exists that maximizes the displacement X.

Figure 4.
Displacement over one cycle X as a function of η under the square gait with stroke amplitudes of (a) ε=π/4 and (b) ε=π/2. Solid curves, numerical integration. Dashed curves, O(ε2) approximation. Dash-dotted ...

We now obtain approximations of the optimal value of η based on the leading-order O(ε2) terms in (3.14). Interestingly, it can be seen from (3.15) that the leading-order expressions for both square and circular gaits are identical up to multiplication by a constant. An important observation is that this relation is fairly general. That is, for any small-amplitude gait and any swimmer's dynamics with the same structure (i.e. not necessarily assuming resistive force theory), the leading-order term of X can always be decomposed into a function of swimmer's geometry multiplied by a function of the gait's unscaled trajectory. This key statement is summarized in the following theorem.

Theorem 3.2 —

Consider a swimmer with two shape variables whose dynamics can be written in the form of equation (2.13), under a small-amplitude gait ϕ(t)=εs(t). The leading-order approximation of the swimmer's displacement X over one period can be written in the form: X(2)=C(s)[center dot]f(η), where C(s) depends on the shape of the unscaled trajectory and f(η) depends on the swimmer's geometric structure.

The theorem, whose proof is detailed in the electronic supplementary material, implies that based on leading order approximation, the optimal geometry is independent of the gait's shape. For our particular swimmer model, the geometry-dependent function in (3.14) is ηl(1−η)3(η+3). As expected, for the cases of η=0 and η=1, the net displacement vanishes, and the optimal value of η is easily obtained as η=0.4101=0.2649, which results in maximal displacement of X*=0.0859ε2 for the square gait and X*=0.0675ε2 for the circular gait. The value of η* is in agreement with the optimal geometry obtained in [21], but not with the one in [13]. This disagreement, which originates from differences in definitions and scaling, is discussed in the sequel.

Next, we use both O(ε2) and O(ε4) terms in (3.14) in order to obtain both optimal ratio η and amplitude ε for the two gaits. The polynomials f2(η) and f4(η) in (3.15) are the same for both gaits up to multiplication by a scalar. However, this is not true for all gaits. For example, it can be shown that for an elliptical gait trajectory, f4(η) is a different polynomial. For a given value of η, (3.14) implies that the locally optimal amplitude ε* for maximal displacement is given by

ε=f2(η)2f4(η).
3.16

Substituting into (3.14) one obtains the optimal displacement as

X=f2(η)4f4(η).
3.17

This implies the existence of an optimal combination of link ratio η and amplitude ε that achieve a local maximum of the displacement. Differentiating (3.17) with respect to η, the optimal value is obtained by finding the roots of a ninth order polynomial in η. As the expression for X* in (3.17) is the same for both square and circular gaits up to a multiplicative constant, the optimal geometry for both gaits is obtained as η*=0.1784. Substituting into (3.16), optimal amplitudes ε and displacements X are obtained for each gait and are presented in table 1. Numerically calculating the optimal geometry η and stroke amplitude ε using Matlab's function fmincon gives the results also presented in table 1.

Table 1.
Values for maximal displacement.

The locally optimal amplitude ε* based on O(ε4) approximation is quite close to the exact value obtained numerically. On the other hand, the large difference in η* can be explained by the large deviation in X for επ/2, as seen in figure 4b. The figure also shows that the O(ε4) expression for X overestimates the exact value. Finally, it is important to note that all optimization results discussed here are limited to finding local rather than global maxima. In fact, continuing to higher order terms in the expansions of (3.14) reveals a minimum of X<0 for ε≈3 rad with larger absolute value than the first positive maximum, and even higher extremum points for larger non-physical values of ε [30].

(d) Comparison with optimal geometry obtained by Becker et al.

As mentioned above, the work of Becker et al. [13] obtained a different value of optimal geometry η. The explanation for this difference is twofold: first, Becker et al. [13] considered maximization of the mean forward velocity under a constant torque at the joint rather than net displacement. Second, they [13] scaled swimming distance by the side link's length l1 rather than the total length l.

We now show that by adopting the definitions of [13] and using our leading-order expressions in proposition 3.1, one obtains optimal geometry that agrees with [13]. First, we consider the case where a constant torque is applied on the rotating joint in each part of the square gait. The relation between the joint torques and the joint velocities is given in (2.16). For the first quarter of the square gait we have that ϕ˙1=0 and so the second row of (2.16) is reduced to τ2=W22ϕ˙2, with Wij the elements of W(ϕ) defined in (2.19). Assuming a constant torque of τ2(t)=τo at the joint for the first quarter cycle of the square gait, we can derive a leading-order expression for the time it would take to complete the quarter gait and then multiply by four in order to obtain

Tτ=dt=4εεdtdϕ2dϕ2=4εεW22(ϕ)τodϕ2=4εεW22(0)τodϕ2+O(ε3)=112l3(1η2)3ctτoε+O(ε3).
3.18

The mean forward velocity V τ=X/Tτ can now be computed to leading order as

Vτ(1)=X(2)Tτ(1)=ε3η(η+3)l2(η+1)3τoct.
3.19

From (3.19), the optimal geometry that maximizes the leading-order mean velocity V (1) is obtained as η*=0.646, which is fundamentally different from the optimal geometry for maximal displacement X obtained above.

Becker et al. [13] also used a different definition of geometric ratio, denoted here as ηb=l0/l1. That is, the length is normalized by that of the side links. The relation between the two definitions of η is given by

ηb=l0l1=2η1η.
3.20

Leading-order expressions for the time of quarter period under constant torque Tτ and the mean velocity V are given in [13] as (after multiplying Tτ by four and adapting some notation)

Tτ=[163(ηb+1)3l13(ηb+4)(ηb+1)2+3ηb+4+O(ε2)](ctl03τo)
3.21

and

Vτ=[34ηb(2ηb+3)(ηb+1)3(ηb+2)l12+O(ε2)](τoctl02).
3.22

Substituting the transformation (3.20), it can be shown that these expressions agree with our derivations in (3.18), (3.19), and also with our leading order expression for net displacement X=V τTτ as given in (3.14). From (3.22), optimal geometry for maximizing V τ has been obtained in [13] as ηb=0.54, which corresponds to η=0.213. The difference from the optimal value of η=0.646 obtained by maximizing V τ in (3.19) is now due to the fact that (3.19) is scaled by l while (3.22) is scaled by l1.

(e) Period time under constant power

We now compute series expansions for the period time Tp under constant power P0=1 as formulated in (2.28), for both square and circular gaits. These expressions are used in the next section for obtaining approximate expressions of Lighthill's energetic efficiency according to equation (2.29). The period time is expanded into a power series as

Tp=εTp(1)+ε2Tp(2)+ε3Tp(3)+.
3.23

Explicit expressions for the different terms Tp(i) depend on the choice of the gait. For the square and circular gaits, these expressions are summarized in the following proposition.

Proposition 3.3 —

For a symmetric, three linked ‘Purcell swimmer’ performing a square or circular gait with amplitude ε, the period times of one full stroke with constant power of Po=1 exerted by the joints are expanded as follows.

Square gait:

Tp=63ctl3(1η2)3ε+6ctl3(1η)4(η+1)(η3+3η23η+1)12(1η2)3ε3+O(ε5).
3.24

Circular gait:

Tp=E(η3+3η23η1η2(η+3))ctl33η2(1η)3(η+3)ε+O(ε3).
3.25

The function E([center dot]) in (3.25) denotes a complete elliptic integral of the second kind, whose definition is reviewed in the electronic supplementary material. The proof of this proposition is also given in the electronic supplementary material.

In order to validate the expansions of Tp in (3.24) and (3.25), figure 5 plots the O(ε) and O(ε3) approximations of Tp as a function of stroke amplitude ε, compared with the exact value of Tp obtained by numerical integration of (2.28), for both gaits (figure 5a, square, and and55b, circular) under equal links length, η=13. It can be seen that for small amplitudes, the period time Tp indeed scales linearly with ε and that for intermediate values of ε the O(ε3) correction term slightly improves the approximation accuracy (square gait only).

Figure 5.
Period time under constant power as a function of stroke amplitude ε for η=13, comparison between exact and approximate expressions. (a) Square gait and (b) circular gait. (Online version in colour.)

4. Approximate expressions for Lighthill's efficiency

Using the expressions found in (3.24),(3.25) for the period times under constant power, along with previous results of net displacement (3.14), we now calculate several approximations of Lighthill's energetic efficiency ξ using (2.29). Then we obtain optimal values of length ratio η and stroke amplitude ε for maximal energetic efficiency, under both square and circular gaits. Using the expansions for X and Tp, the efficiency ξ can be written as

ξ=X2Tp2=ε2(X(2)+ε2X(4)+)2(Tp(1)+ε2Tp(3)+)2.
4.1

Calculating Taylor expansion of (4.1) in powers of ε, one obtains

ξ=ε2(X(2))2(Tp(1))2+2ε4X(2)Tp(1)X(4)X(2)Tp(3)(Tp(1))3+O(ε6).
4.2

Explicit approximations for the efficiency ξ under square and circular gaits are discussed below.

(a) Square gait

Substituting (3.14), (3.15) and (3.24) into (4.2), the first two elements in the expansion of Lighthill's energetic efficiency under the square gait are given by

ξ(2)=(X(2))2(Tp(1))2=3η2(1η)3(η+3)232ctl(η+1)3
4.3

and

ξ(4)=η2ε4(η1)3(η+3)(η9+5η83η739η637η5+71η4+263η3+371η248η+56)256ctl(η+1)5.
4.4

Figure 6a plots the O(ε2) and O(ε4) approximations of ξ as a function of ε in dashed and dash-dotted curves, respectively, under the square gait with length ratio of η=13. For comparison, the exact value of ξ obtained from numerical integration is plotted in solid curve. While the O(ε2) approximation grows monotonically with ε, adding the next-order correction of O(ε4) also captures a local maximum in the efficiency, which is obtained at an intermediate stroke amplitude of ε≈1. Nevertheless, this approximation of ξ becomes negative for larger values of ε, which is non-physical according to the definition of ξ in (2.29). A more reasonable approximation is obtained by plugging the O(ε4) and O(ε3) approximations for X and Tp, respectively, into (2.29). With a slight abuse of notation, this approximation is denoted here as O(ε4/3), and is given by

ξ(4/3)=(ε2X(2)+ε4X(4))2(εTp(1)+ε3Tp(3))2.
4.5

Note that this approximation is different from the direct expansion of (4.1) into (4.2), as demonstrated below. A plot of this approximation is overlayed as a dotted line in figure 6a. It can be seen that the O(ε4/3) approximation is always positive, and correctly captures the general trend of ξ as a function of ε. Finally, it can be seen from the plot that the exact value of ξ attains a higher maximum value at ε≈3, where the direction of swimming is reversed. This second maximum is not captured by any of the approximate expressions mentioned above. Nevertheless, as already mentioned in [13], this maximum at large strokes where the joint angles approach ±π is often impractical, due to possible collisions between the links. Moreover, resistive force theory which assumes negligible hydrodynamic interaction between the links is no longer valid in this range [30]. Thus, we focus here on approximations of the first maximum of ξ which is attained at moderate amplitudes, as discussed next.

Figure 6.
Lighthill's energetic efficiency ξ under the square gait: (a) plot of ξ versus ε for η=13. (b) Plot of ξ versus η for ε=1. Solid curves, exact (numeric) computation. Dashed, dash-dotted and dotted ...

(i) Efficiency optimization (square gait)

We now use the approximations of Lighthill's energetic efficiency in order to obtain locally optimal values of both length ratio η and stroke amplitude ε. Using only the leading-order approximation of the efficiency in (4.3) gives an optimal length ratio of η*=0.327. Taking also the next-order correction term into account and implementing the same calculation steps shown above in §3c, the optimal geometry and amplitude are obtained as

η=0.2853,ε=0.908rad=52andξ=0.5838.
4.6

Using O(ε(4/3)) approximation, the optimal geometry and amplitude are obtained by solving a system of two polynomial equations as

η=0.2754,ε=1.1038rad=63.2andξ=0.7294.
4.7

Finally, numerical calculation of ξ under the square gait and conducting optimization using Matlab's function fmincon gives optimal values of

η=0.2879,ε=1.133rad=64.9andξ=0.7709.
4.8

Figure 6b shows a plot of the approximations of ξ as a function of η for a large amplitude of ε=1, compared with the exact computation of ξ obtained numerically. While there are large discrepancies in the value of the efficiency ξ (best captured by the O(ε4/3) approximation), all approximations predict optimal values around η=0.28.

(b) Circular gait

For the circular gait, the constant-power period time Tp has been approximated only to first order O(ε) in (3.25). Using (3.14), the expansion of Lighthill's efficiency in (4.2) can thus be obtained only to leading order, as

ξ=ε23(1η)3(η+3)ctlE(η3+3η23η1η2(η+3))+O(ε3).
4.9

Similar to (4.5), we define the O(ε4/1) approximation as

ξ(4/1)=(ε2X(2)+ε4X(4))2(Tp(1))2.
4.10

Substituting (3.14) and (3.25) into (4.10) gives the expression for the O(ε4/1) approximation (not shown for brevity). The efficiency ξ and its approximations are plotted in figure 7a as a function of the stroke amplitude ε under equal links lengths η=13, and in figure 7b as a function of η under amplitude of ε=1. In both plots, O(ε2) approximations appear in dashed curves, O(ε4/1) approximations appear in dotted curves, and the exact values computed numerically are shown in solid curves. As before, the leading-order O(ε2) approximation does not capture the optimum with respect to amplitude but both approximations capture optimum with respect to η.

Figure 7.
Lighthill's energetic efficiency ξ under the circular gait: (a) plot of ξ versus ε for η=13. (b) Plot of ξ versus η for ε=1. Solid curves, exact (numeric) computation. Dashed and dotted curves, approximations ...

(ii) Efficiency optimization (circular gait)

Using the leading-order approximation for ξ in (4.9) and numerically searching for maximum with respect to η, optimal geometry is obtained as η*=0.3139, with an efficiency of ξ(2)=1.5565. Using the O(ε4/1) approximation in (4.10), optimal values of both ε and η are obtained numerically as

η=0.2379,ε=1.3821radandξ=1.22.
4.11

Finally, numerical calculation of ξ under the circular gait and conducting optimization using Matlab's fmincon function give optimal values of

η=0.2684,ε=1.3747radandξ=1.2777.
4.12

Again, the O(ε4/1) expression in (4.10) achieves a reasonable approximation of the optimum, as also seen from figure 7a,b.

An important observation is that even in leading-order O(ε2) approximation of Lighthill's energetic efficiency ξ, the optimal geometry η* for maximizing ξ actually depends on the shape of the chosen gait trajectory. This can clearly be seen in the different dependence of ξ(2) on η in (4.3) and (4.9). This is a substantial difference from the case of maximizing displacement X, where dependence of the leading-order approximation in η is independent of the gait's shape, as manifested in theorem 3.2.

5. Conclusion

In this work, we have analysed the motion of the three-link ‘Purcell swimmer’. We provided a systematic method for deriving an expansion of the velocities of the swimmer and presented leading-order expressions and next-order corrections for the net displacement over one period in the cases of a square and circular gait. Examination of the correction terms confirms that there is a reversal in the direction of the net displacement at high amplitudes, a result which has previously only been shown numerically. The gait amplitude and swimmer geometry that optimize the displacement over one period were approximated using the first two terms in the expansion. Additionally, by writing asymptotic expressions for the period time under constant power expenditure, we were able to write, for the first time, leading-order term for the energetic efficiency of the square and circular gaits as well as a next-order correction term for the square gait. Once again, we used the obtained expressions in order to find the energetically optimal gait amplitude and swimmer geometry. The results demonstrate the utility of perturbation methods for obtaining approximate explicit expressions for nonlinear dynamics of locomotion systems, which enable analysis and optimization of their performance.

We now briefly discuss some limitations of this work and list possible directions for future extension of the research. First, the swimmer's dynamic equations have been formulated using the simplification of resistive force theory [26,27], in which hydrodynamic interaction between the links is neglected. It is well known (cf. [15,31]) that this assumption holds only for highly slender links, and for small stroke amplitudes where the gap between the links remains large even in the vicinity of the joints. Hydrodynamic interactions may be accounted for by using more refined models of slender body theory as in [30,18]. The decoupled relations in (2.5) should then be modified to account for inter-link resistance, and the analysis will probably become purely numerical due to the added complexity of the dynamics. Nevertheless, the structure and geometric symmetries of the dynamics in (2.11) will be maintained without qualitative changes. In fact, according to the numerical investigation in [30,18], no large quantitative changes in the results are expected even for moderate stroke amplitudes (ε≈2 rad).

Second, the gait optimization conducted in this work was limited to varying the amplitude of predefined shape trajectories and did not address the possibility of other trajectories which may have better performance. Shape optimization of swimmers has been extensively studied, cf. [3234]. One possible approach [18,35,36] is to represent a subset of all possible periodic shape changes by a finite set of variables (e.g. coefficients of a truncated Fourier series) and perform numerical optimization over this discrete set of variables. Another approach is optimal control theory [37], which is based on calculus of variations. This approach has been exploited by Alouges et al. [38] for obtaining energy-optimal gaits of unidirectional axisymmetric swimmers and by Giraldi et al. [21] for gait optimization of multi-link microswimmers under input constraints. In our recent work [39], we have used optimal control in order to obtain unconstrained optimal gaits for maximal displacement of Purcell's swimmer, which exactly reproduce the gaits obtained numerically in [18]. Furthermore, we are currently working on extending the results of [38] to planar motion of multi-link microswimmers.

Third, experiments conducted on a macro-scale swimmer prototype in a highly viscous fluid [14] have recently demonstrated the symmetry properties of Purcell's swimmer performing certain symmetric gaits. Similar experiments can be done to verify some results from this paper, e.g. the dependency of net displacement on swimmer geometry and gait amplitude. Finally, in many practical situations, the actuation of robotic swimmers will likely be by applying torques at the joints rather than prescribing the joint angles directly [40]. This calls for the formulation of an asymptotic expansion of the relation between the joint torques and the joint angles.

Supplementary Material

Supplementary document:

Authors' contributions

Both authors contributed equally to the writing of this paper and both authors gave final approval for publication.

Competing interests

We declare we have no competing interests.

Funding

This work is supported by the Israel Science Foundation (ISF) under grant no. 567/14.

References

1. Lauga E, Powers TR 2009. The hydrodynamics of swimming microorganisms. Rep. Prog. Phys. 72, 096601 (doi:10.1088/0034-4885/72/9/096601)
2. Lauga E.2016. Bacterial hydrodynamics. Annu. Rev. Fluid Mech. (http://arxiv.org/abs/1509.02184. )
3. Abbott JJ, Lagomarsino MC, Zhang L, Dong L, Nelson BJ 2009. How should microrobots swim? Int. J. Robot. Res. 28, 1434–1447. (doi:10.1177/0278364909341658)
4. Kósa G, Jakab P, Hata N, Jólesz F, Neubach Z, Shoham M, Zaaroor M, Székely G 2008. Flagellar swimming for medical micro robots: theory, experiments and application. In 2nd IEEE RAS & EMBS Int. Conf. on Biomedical Robotics and Biomechatronics (BioRob), pp. 258–263. New York, NY: IEEE.
5. Jang B. et al. 2015. Undulatory locomotion of magnetic multilink nanoswimmers. Nano Lett. 15, 4829–4833. (doi:10.1021/acs.nanolett.5b01981) [PubMed]
6. Gao W. et al. 2012. Cargo-towing fuel-free magnetic nanoswimmers for targeted drug delivery. Small 8, 460–467. (doi:10.1002/smll.201101909) [PubMed]
7. Peyer KE, Zhang L, Nelson BJ 2013. Bio-inspired magnetic swimming microrobots for biomedical applications. Nanoscale 5, 1259–1272. (doi:10.1039/C2NR32554C) [PubMed]
8. Sitti M, Ceylan H, Hu W, Giltinan J, Turan M, Yim S, Diller E 2015. Biomedical applications of untethered mobile milli/microrobots. Proc. IEEE 103, 205–224. (doi:10.1109/JPROC.2014.2385105)
9. Happel J, Brenner H 1965. Low Reynolds number hydrodynamics. Englewood Cliffs, NJ: Prentice-Hall.
10. Purcell EM. 1977. Life at low Reynolds number. Am. J. Phys. 45, 3–11. (doi:10.1119/1.10903)
11. Berg HC, Anderson RA 1973. Bacteria swim by rotating their flagellar filaments. Nature 245, 380 (doi:10.1038/245380a0) [PubMed]
12. Taylor G. 1951. Analysis of the swimming of microscopic organisms. Proc. R. Soc. Lond. A 209, 447–461. (doi:10.1098/rspa.1951.0218)
13. Becker LE, Koehler SA, Stone HA 2003. On self-propulsion of micro-machines at low Reynolds number: Purcell's three-link swimmer. J. Fluid Mech. 490, 15–35. (doi:10.1017/S0022112003005184)
14. Gutman E, Or Y 2016. Symmetries and gaits for Purcell's three-link microswimmer model. IEEE Trans. Robot. 32, 53–69. (doi:10.1109/TRO.2015.2500442)
15. Avron JE, Raz O 2008. A geometric theory of swimming: Purcell's swimmer and its symmetrized cousin. New J. Phys. 10, 063016 (doi:10.1088/1367-2630/10/6/063016)
16. Hatton RL, Ding Y, Choset H, Goldman DI 2013. Geometric visualization of self-propulsion in a complex medium. Phys. Rev. Lett. 110, 078101 (doi:10.1103/PhysRevLett.110.078101) [PubMed]
17. Hatton RL, Choset H 2013. Geometric swimming at low and high Reynolds numbers. IEEE Trans. Robot. 29, 615–624. (doi:10.1109/TRO.2013.2251211)
18. Tam D, Hosoi AE 2007. Optimal stroke patterns for Purcell's three-link swimmer. Phys. Rev. Lett. 98, 068105 (doi:10.1103/PhysRevLett.98.068105) [PubMed]
19. Shapere A, Wilczek F 1989. Geometry of self-propulsion at low Reynolds number. J. Fluid Mech. 198, 557–585. (doi:10.1017/S002211208900025X)
20. Nayfeh AH. 2008. Perturbation methods. New York, NY: John Wiley & Sons.
21. Giraldi L, Martinon P, Zoppello M 2015. Optimal design of Purcell's three-link swimmer. Phys. Rev. E 91, 023012 (doi:10.1103/PhysRevE.91.023012) [PubMed]
22. Passov(Gutman) E, Or Y 2012. Dynamics of Purcell's three-link microswimmer with a passive elastic tail. Eur. Phys. J. E 35, 1–9. (doi:10.1140/epje/i2012-12001-6) [PubMed]
23. Gutman E, Or Y 2014. Simple model of a planar undulating magnetic microswimmer. Phys. Rev. E 90, 013012 (doi:10.1103/PhysRevE.90.013012) [PubMed]
24. Gutman E, Or Y 2016. Optimizing an undulating magnetic microswimmer for cargo towing. Phys. Rev. E 93, 063105 (doi:10.1103/PhysRevE.93.063105) [PubMed]
25. Man Y, Lauga E 2013. The wobbling-to-swimming transition of rotated helices. Phys. Fluids 25, 071904 (doi:10.1063/1.4812637)
26. Cox RG. 1970. The motion of long slender bodies in a viscous fluid part 1. General theory. J. Fluid Mech. 44, 791–810. (doi:10.1017/S002211207000215X)
27. Gray J, Hancock G 1955. The propulsion of sea-urchin spermatozoa. J. Exp. Biol. 32, 802–814.
28. Childress S. 2012. A thermodynamic efficiency for Stokesian swimming. J. Fluid Mech. 705, 77–97. (doi:10.1017/jfm.2011.561)
29. Lighthill J. 1975. Mathematical biofluiddynamics. Philadelphia, PA: SIAM.
30. Huber G, Koehler SA, Yang J 2011. Micro-swimmers with hydrodynamic interactions. Math. Comput. Model. 53, 1518–1526. (doi:10.1016/j.mcm.2010.04.002)
31. Alouges F, DeSimone A, Heltai L, Lefebvre-Lepot A, Merlet B 2013. Optimally swimming stokesian robots. Discrete Continuous Dyn. Syst. B 18, 1189–1215. (doi:10.3934/dcdsb.2013.18.1189)
32. Delfour MC, Zolésio J-P 2011. Shapes and geometries: analysis, differential calculus, and optimization, 2nd edn Advances in design and control, vol. 4 Philadelphia, PA: SIAM.
33. Moubachir M, Zolésio J-P 2006. Moving shape analysis and control: applications to fluid structure interactions. Pure and applied mathematics, vol. 277 London, UK: Chapman and Hall.
34. Walker SW. 2015. The shapes of things: a practical guide to differential geometry and the shape derivative. Advances in design and control, vol. 28 Philadelphia, PA: SIAM
35. Avron JE, Gat O, Kenneth O 2004. Optimal swimming at low Reynolds numbers. Phys. Rev. Lett. 93, 186001 (doi:10.1103/PhysRevLett.93.186001) [PubMed]
36. Shapere A, Wilczek F 1989. Efficiencies of self-propulsion at low Reynolds number. J. Fluid Mech. 198, 587–599. (doi:10.1017/S0022112089000261)
37. Bryson A, Ho Y-C 1975. Applied optimal control: optimization, estimation, and control. Washington DC: Hemisphere Pub. Corp.
38. Alouges F, DeSimone A, Lefebvre A 2009. Optimal strokes for axisymmetric microswimmers. Eur. Phys. J. E 28, 279–284. (doi:10.1140/epje/i2008-10406-4) [PubMed]
39. Wiezel O, Or Y In preparation Using optimal control to obtain maximum displacement gait for Purcell's three-link swimmer. In Proc. 2016 IEEE Conf. on Decision and Control (CDC).
40. Or Y. 2012. Asymmetry and stability of shape kinematics in microswimmers’ motion. Phys. Rev. Lett. 108, 258101 (doi:10.1103/PhysRevLett.108.258101) [PubMed]

Articles from Proceedings. Mathematical, Physical, and Engineering Sciences are provided here courtesy of The Royal Society