## 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.

## 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 [3–5]. Such mini-robots may have many applications in medicine, performing medical procedures in a minimally invasive way and delivering drugs with high precision [6–8]. 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*≪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 2*a*). 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 2*b*), at which the two joint angles oscillate sinusoidally with a quarter-period phase shift.

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 2*b*); (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 *l*_{0},*l*_{1},*l*_{2}, with *l*_{1}=*l*_{2}, *l*=*l*_{0}+*l*_{1}+*l*_{2} and *η*=*l*_{0}/*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

*x*−

*y*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

*i*th link is described by the linear velocity of its centre

**v**

_{i}and the link's angular velocity

*ω*

_{i}, which are augmented in the vector

**V**

_{b}=(

**v**

_{b},

*ω*

_{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

*I*

_{ij}, such that

*I*

_{ij}=±1 if the location of the

*i*th link with respect to the body frame is affected by the

*j*th joint, with the sign determined by the direction of the joint axis.

*I*

_{ij}=0, if the

*j*th joint does not affect the

*i*th 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
**r**_{i} is the position of the centre of the *i*th link and **b**_{j} is the position of the *j*th joint. In matrix form, the velocity **V**_{i} is related to the body velocity **V**_{b} and shape velocity **T**_{i}(**q**,** ϕ**) and

**E**

_{i}(

**q**,

**) for**

*ϕ**i*=0,1,2 are given by

*α*

_{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

**f**

_{i}and torque

*m*

_{i}on the

*i*th 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:

*i*th link and

*μ*is the dynamic viscosity of the fluid,

*a*is the radius of the links and

*l*

_{c}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

*i*th link as

**F**

_{i}=(

**f**

_{i},

*m*

_{i}), the relation (2.4) can be written in matrix form

**T**

_{i}from the kinematic relation (2.2) and augmenting in a vector

**F**

_{b}=(

**f**

_{b},

*m*

_{b}), (2.7) is written in matrix form as

**q**implies that the body velocity satisfies

**F**

_{b}=0. Substituting this into (2.10) then yields the nonlinear differential equations that govern the swimmer's dynamics:

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])
*θ*. The matrix **G**(** ϕ**) obeys some symmetry relations due to the swimmer's structure (see [14] for details):

*ϕ*

_{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

*j*th joint, these torques are balanced by hydrodynamic forces and torques

**f**

_{i}and

*m*

_{i}, giving rise to the following relation:

**E**

_{i}:

**=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*)=

*εs*

_{i}(

*t*), where

*s*

_{i}represents the ‘unscaled’ shape trajectory. Let us also define the vector

**s**=(

*s*

_{1},

*s*

_{2})

^{T}. For the cases of square and circular gaits, the joint angles can be written in the unscaled form as

*ε*by setting

*ϕ*

_{i}(

*t*)=

*εs*

_{i}(

*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

*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

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 ** ϕ**(

*σ*(

*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

*σ*(

*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

The period time *T* under a constant mechanical power *P*=*P*_{o} can be found using the following derivation. Consider the gait ** ϕ**(

*σ*)=

**(**

*ϕ**σ*(

*t*)), where

*σ*∈[

*σ*

_{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

*F*(

*σ*)=

**′(**

*ϕ**σ*)

**W**(

*σ*)

**′(**

*ϕ**σ*) we have

**′(**

*ϕ**σ*)≠0 and d

*σ*/d

*t*>0, for all

*σ*and

*t*. The period time of a gait under constant mechanical power

*P*=

*P*

_{o}is thus obtained as

*P*

_{o}cancels out and so the calculations can be done for

*P*

_{o}=1 without loss of generality. Thus, Lighthill's efficiency of a gait

**(**

*ϕ**σ*) is uniquely given as

*T*

_{p}now denotes the period time under constant power of

*P*

_{o}=1, as given in (2.28), and the constants

*c*

_{t}

*l*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
*θ*, and thus the ODE for *θ* can be written as
*θ* and *ϕ*_{i}(*t*)=*εs*_{i}(*t*):
*ε*.
**G**_{31}, **G**_{32} are even functions in (*ϕ*_{1},*ϕ*_{2}) and so *θ*^{(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
**G** in (3.6). Expanding equation (2.13) by substituting expansions for *θ*(*t*), **D**(*θ*), **G**(** ϕ**) and also of

**q**and

*ϕ*

_{i}(

*t*)=

*εs*

_{i}(

*t*) and collecting orders of

*ε*, one obtains

**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*

*With*

*where for the square gait given in*(2.21),

*we have*

*and for the circular gait given in*(2.22)

*C*

_{2}=

*π*/16,

*C*

_{4}=

*π*/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 *a*, square, and 3*b*, 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}).

### (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 4*a*,*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*.

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***)⋅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 *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 *f*_{2}(*η*) and *f*_{4}(*η*) 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, *f*_{4}(*η*) is a different polynomial. For a given value of *η*, (3.14) implies that the locally optimal amplitude *ε** for maximal displacement is given by
*η* 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.

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 4*b*. 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 *l*_{1} 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 *W*_{ij} 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

*V*

_{τ}=

*X*/

*T*

_{τ}can now be computed to leading order as

*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}=*l*_{0}/*l*_{1}. That is, the length is normalized by that of the side links. The relation between the two definitions of *η* is given by
*T*_{τ} and the mean velocity *V* are given in [13] as (after multiplying *T*_{τ} by four and adapting some notation)
*X*=*V* _{τ}*T*_{τ} as given in (3.14). From (3.22), optimal geometry for maximizing *V* _{τ} has been obtained in [13] as *η*=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 *l*_{1}.

### (e) Period time under constant power

We now compute series expansions for the period time *T*_{p} under constant power *P*_{0}=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

### 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* *P*_{o}=1 *exerted by the joints are expanded as follows*.

*Square gait*:
*Circular gait*:

The function *E*(⋅) 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 *T*_{p} in (3.24) and (3.25), figure 5 plots the *O*(*ε*) and *O*(*ε*^{3}) approximations of *T*_{p} as a function of stroke amplitude *ε*, compared with the exact value of *T*_{p} obtained by numerical integration of (2.28), for both gaits (figure 5*a*, square, and 5*b*, circular) under equal links length, *T*_{p} indeed scales linearly with *ε* and that for intermediate values of *ε* the *O*(*ε*^{3}) correction term slightly improves the approximation accuracy (square gait only).

## 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 *T*_{p}, the efficiency *ξ* can be written as
*ε*, one obtains
*ξ* 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

Figure 6*a* 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 *ξ* 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 *T*_{p}, respectively, into (2.29). With a slight abuse of notation, this approximation is denoted here as *O*(*ε*^{4/3}), and is given by
*a*. 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.

#### (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
*O*(*ε*^{(4/3)}) approximation, the optimal geometry and amplitude are obtained by solving a system of two polynomial equations as
*ξ* under the square gait and conducting optimization using Matlab's function **fmincon** gives optimal values of

Figure 6*b* 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 *T*_{p} 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
*O*(*ε*^{4/1}) approximation as
*O*(*ε*^{4/1}) approximation (not shown for brevity). The efficiency *ξ* and its approximations are plotted in figure 7*a* as a function of the stroke amplitude *ε* under equal links lengths *b* 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 *η*.

#### (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
*ξ* under the circular gait and conducting optimization using Matlab's **fmincon** function give optimal values of
*O*(*ε*^{4/1}) expression in (4.10) achieves a reasonable approximation of the optimum, as also seen from figure 7*a*,*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. [32–34]. 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.

## 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.

- Received June 2, 2016.
- Accepted July 15, 2016.

- © 2016 The Author(s)

Published by the Royal Society. All rights reserved.