## Abstract

We describe a boundary-element method used to model the hydrodynamics of a bacterium propelled by a single helical flagellum. Using this model, we optimize the power efficiency of swimming with respect to cell body and flagellum geometrical parameters, and find that optima for swimming in unbounded fluid and near a no-slip plane boundary are nearly indistinguishable. We also consider the novel optimization objective of torque efficiency and find a very different optimal shape. Excluding effects such as Brownian motion and electrostatic interactions, it is demonstrated that hydrodynamic forces may trap the bacterium in a stable, circular orbit near the boundary, leading to the empirically observable surface accumulation of bacteria. Furthermore, the details and even the existence of this stable orbit depend on geometrical parameters of the bacterium, as described in this article. These results shed some light on the phenomenon of surface accumulation of micro-organisms and offer hydrodynamic explanations as to why some bacteria may accumulate more readily than others based on morphology.

## 1. Introduction

Monotrichous bacteria, such as *Rhodobacter sphaeroides* and *Vibrio alginolyticus*, have a single flagellum at one of the poles of the cell body. They are well suited to mathematical modelling due to geometrical simplicity. The flagellum is approximately helical and extends directly behind the cell body during periods of steady swimming. Movies of *R. sphaeroides* (available from Howard Berg’s laboratory website, http://www.rowland.harvard.edu/labs/bacteria/index_movies.html) indicate that, to a good approximation, the flagellum appears rigid and simply rotates as the organism swims. However, flagella can exhibit elastic properties, for example, undergoing conformational changes during the reorientation phase of bacterial motility (Armitage *et al.* 1999).

When examining the fluid dynamics of bacterial propulsion, it is found that the Reynolds number is very low (*Re*≈10^{−5}) and so fluid flow may be described by the incompressible Stokes equations. Linearity of these equations allows the superposition of singular solutions to construct flows that satisfy the required boundary conditions, a technique widely used in studies to date. However, the microscopic scale at which bacteria live also means that Brownian motion is potentially an issue. Fortunately, the results of Dusenbery (1997) suggest that neglecting Brownian motion is justified for bacteria as large as *R. sphaeroides*, with a length of 1.73 μm (Slovak *et al.* 2005). This is also evident from inspection of the bacterial swimming movies produced by Berg’s laboratory.

The movement of such flagellated micro-organisms through bulk liquid is well understood, both on a statistical level (where, for example, run-and-tumble trajectories can be thought of as random walks) and on the individual level. Mathematical modelling of swimming in low Reynolds numbers began with Taylor (1951, 1952) and centred on slender-body theory (SBT; Burgers 1938; Hancock 1953), as well as the closely related resistive-force theory (RFT), first used by Gray & Hancock (1955) in their calculations for planar bending waves propagating along flagella of sea-urchin spermatozoa. These theories exploit the slenderness of flagella to reduce the problem to one dimension. However, accuracy is compromised when the flagellum cannot reasonably be treated as a series of locally straight, isolated filament sections in an infinite body of fluid, for example, in near-surface swimming or if the flagellum is tightly wound. Many analytical and numerical results have been presented for spherical cell bodies and helical or sinusoidal, planar flagellar beats in unbounded fluids (Chwang & Wu 1971; Keller & Rubinow 1976; Lighthill 1976; Higdon 1979*a*,*b*).

While RFT is a simple local approximation to SBT, it agrees very well with SBT analysis and adequately predicts the dynamics of swimming in unbounded fluid, except where the cell body is large and has significant influence (Johnson & Brokaw 1979). However, even SBT only approximately satisfies the boundary conditions of the model and there are more accurate, though more sophisticated, methods. Phan-Thien *et al.* (1987) used a boundary-element method (BEM) to analyse the motion of an ellipsoidal cell body propelled by a rigidly rotating helical flagellum. Good agreement was found when compared to the SBT results of Higdon (1979*b*) for swimming in unbounded fluids. More recently, Fujita & Kawai (2001) have drawn the same conclusions with quantitatively very similar results in their own BEM study.

Authors, such as Higdon (1979*b*), Phan-Thien *et al.* (1987) and Fujita & Kawai (2001), have addressed the question of optimizing swimming speed in bulk fluid with respect to shape using power efficiency as the objective function. While this may have some relevance, particularly in the design of artificial microswimmers, there are a host of other factors that would affect biological fitness of bacteria in nature (Young 2006). Even if the only aim was to make the fastest possible swimmer, there would be several conditions to consider, such as the flagellar motor capabilities, as well as the cell and flagellum shapes.

In bacterial propulsion, the driving force comes from a rotary motor embedded in the cell membrane and connected to the flagellar filament by a short, flexible segment called the hook. The motor is normally situated at a pole of the body on monotrichous bacteria. Studies have found that the torque produced by the motor is relatively constant over a large range of rotational frequencies before it drops roughly linearly to zero at higher frequencies (Berg & Turner 1993; Chen & Berg 2000; Xing *et al.* 2006). The low-frequency torque, the frequency at which the torque begins to decrease and the slope of decrease all vary between species, and some properties depend on environmental conditions such as temperature (Sowa & Berry 2008).

The motor clearly has an important role in determining swimming speeds and the fastest swimmer powered by a given motor is not necessarily the most power efficient. Given that the torque produced by motors is approximately constant at lower rotation rates, we will examine torque efficiency as an alternative optimization criterion. This finds the shape of the fastest swimmer for a given motor torque, and may be relevant when the swimmer is limited by motor output rather than available energy. It is known that motor frequencies are within the approximately constant torque regime during swimming in some species, though not all (Li & Tang 2006). Frequencies may also enter or leave this regime depending on environmental conditions such as temperature and viscosity (Lowe *et al.* 1987; Sowa *et al.* 2003). Our first aim will be to contrast the result of torque optimization with the power optimum and with examples from observation, allowing an exploration of whether optimization studies are robust to the influence of additional biology in the choice of objective functions.

There is a growing interest in developing our understanding of how micro-organisms swim through confined environments, such as in narrow channels or porous media (Biondi *et al.* 1998; Galajda *et al.* 2007). This is relevant in a wide range of scenarios. Bacteria moving through a host organism or in microdevices inevitably experience the effects of solid boundaries, which may even be exploited to achieve certain aims, such as sorting cells by size (Hulme *et al.* 2008). Behaviour near surfaces is also important for initiating biofilm formation (O’Toole *et al.* 2000) and modulating bacterial transport in bioremediation of contaminated groundwater (Sherwood *et al.* 2003). It also affects the foraging patterns of organisms that feed on nutrients diffusing from surfaces. In such cases, it would be beneficial to remain close to the surface so these organisms may exhibit features that encourage this. As an example, the marine bacterium *V. alginolyticus* adopts a run-and-reverse swimming strategy rather than run-and-tumble when it is near a large source of food. It has been found that backward swimming increases their residence time near solid surfaces. On the other hand, tumbling could result in the bacterium swimming away from the surface, which is potentially a source of nutrients (Magariyama *et al.* 2005).

Swimming near no-slip plane boundaries was analysed as early as Reynolds (1965), where the swimmer was an infinite waving sheet. RFT and SBT have been extended to half-space domains (Blake 1974; Katz *et al.* 1975) and there have also been SBT studies on flagellar motility near surfaces (Smith *et al.* 2009). Separations from the wall of the order of the slender-body length can be accommodated (Barta & Liron 1988), but the fundamental assumptions of SBT mean that the model has not been formally justified when a freely swimming organism is much closer than this. Nevertheless, post hoc testing typically shows that the numerical residuals are relatively small, provided the wall separation is greater than a few flagellum radii.

BEMs can model close surfaces more accurately and were used by Ramia *et al.* (1993) to give a comparison of the kinematics of swimmers in unbounded fluid, near a plane boundary, between two parallel-plane boundaries and near another swimmer. The path followed by a swimmer near a single-plane boundary was found to be circular, as is commonly observed under the microscope. This effect has been explained (Lauga *et al.* 2006) and reported in a number of numerical studies. However, any long-time-scale simulations have invariably had to terminate as the swimmer approached and descended into the boundary (as reported by Ramia *et al.* 1993, for example). Apart from curving the trajectory, the boundary also seems to attract the swimmer (Berke *et al.* 2008).

Short-ranged forces between the swimmer and the boundary, such as van der Waals, electrostatic and steric interactions, have been hypothesized to play an important role in trapping microswimmers near surfaces (Frymier *et al.* 1995). Extra interaction potentials are necessary in simulations to account for surface effects at these close distances, creating a finite preferred separation from the boundary (Li *et al.* 2008). However, two different kinds of ‘trapping’ at surfaces have been documented (Vigeant *et al.* 2002): cells can swim freely but at a fixed distance from the surface for extended periods of time or they can actually adhere to the surface and become immobilized. Evidence suggests that short-ranged forces are only responsible for the latter kind of entrapment, whereas hydrodynamic interactions alone are sufficient to draw cells to swim along surfaces. We will focus on the mobile type of entrapment, which Vigeant *et al.* (2002) suggest occurs beyond the range of the additional surface forces.

As mentioned above, attraction of swimmers towards the boundary has already been explained by hydrodynamics. Stable trajectories at a finite separation from a boundary appear to be more elusive, but have been shown to exist purely as a result of hydrodynamics. In the study by Smith *et al.* (2009), the swimmers were spermatozoa and both planar and ‘elliptical helicoid’ beat patterns gave rise to stable boundary accumulation.

Goto *et al.* (2005) deduced the possibility of motion at constant separation from a wall for bacteria by showing that pitching motion is stable when swimming forwards roughly parallel to the boundary. However, they did not seek a stable combination of height and pitch angle, nor did they perform numerical tracking of a swimmer to demonstrate stable motion. Thus, our second aim will be to verify this surface accumulation for monotrichous bacteria and characterize the geometrical aspects of swimmers that encourage or prevent stable surface motility.

We will develop a boundary-element model for a simple bacterial cell, assuming that the surface of the body and flagellum is smooth (i.e. free of pili and other irregularities). As motivated above in the context of *R. sphaeroides*, we neglect Brownian motion, short-range surface forces and structural deformations. The geometry of our bacterial model and the dynamics we consider will be explained further in §2, along with some techniques for analysing the resulting behaviour. Our simulation results will be described in §3 and some interpretations and implications of our findings in the context of questions of bacterial boundary accumulation and optimal bacterial swimming will be given in §4. We will then summarize the investigation in §5.

## 2. Methods

### (a) The boundary-integral equation

BEMs are well established for numerical investigations of a range of physical phenomena, including flagellar motility in viscous fluids (Phan-Thien *et al.* 1987; Goto *et al.* 2001; Smith *et al.* 2009). The present work follows the approach in Pozrikidis (2002).

The unforced, incompressible Stokes equations are
2.1
where *p* and ** u** are the pressure and flow-velocity fields, respectively, and

*μ*is the dynamic viscosity of the fluid.

The flow in a force-free fluid domain, *V* , due to a rigidly moving object or collection of such objects with smooth boundary ∂*V* can be expressed as
2.2
where ** G** is the Stokeslet Green function in the appropriate domain and

**is the traction distribution over the object’s boundary. Expressions for the Stokeslet are given in §1 of the electronic supplementary material. From the boundary-integral equation (2.2), we form a system of linear equations relating the velocity at discrete collocation points on the boundary to the forces at those points. Gauss–Legendre quadrature is used to integrate quantities over curved triangular surface elements following suggestions from Pozrikidis (2002).**

*f*### (b) Model geometry

The model for the bacterium consists of two rigid bodies: an ellipsoid representing the cell body and a thin cylinder with hemispherical ends, curved into a helical shape with an amplitude envelope suggested by Higdon (1979*b*) so that the point of attachment is on the axis of the helix. The centre line of the flagellum is given by
2.3
, where is the extent of the flagellum along the tail axis. The helical amplitude and wavenumber are *a* and *k*, respectively, and *k*_{E} is the amplitude envelope factor. The helical amplitude grows from zero to its full value over a region of length roughly 2/*k*_{E}, which we refer to as the flagellum starting region. Please refer to figure 1 and table 1 for an illustration of the model bacterium and an explanation of geometrical parameters used throughout. This model is similar to that used by Phan-Thien *et al.* (1987), but with a refinement of the filament ends. Isoparametric, quadratic triangular elements are used for greater accuracy in approximating the curved surfaces and surface distributions. Details of the bacterial surface mesh can be found in §2 of the electronic supplementary material.

For all simulations reported here, the cell body–flagellum junction is at a pole of the body (as depicted in figure 1). This is the normal configuration for monotrichous bacteria. In the model, there is a small separation (of the order of a flagellum radius) between the cell body and the flagellum to reduce numerical errors that may arise when two nearby surfaces have different velocities. The reference point of the bacterium is taken to be the starting tip of the flagellum, which we call the junction position, *x*^{J}. Since the relative motion of the flagellum with respect to the cell body is confined to rotations about its axis, the entire configuration may be described in the stationary reference frame by the junction position, *x*^{J}, the head orientation vector, *e*^{H}, and the phase of the tail about its axis, *ϕ*^{T}.

### (c) Solving instantaneous dynamics

This model supposes that the bacterium’s cell body has a certain translational and rotational velocity, ** U** and

*Ω*^{H}, respectively, with the rotational velocity taken about the point

*x*^{J}. The flagellum additionally rotates relative to the body at a rate

*Ω*

^{T}≡

*d*

*ϕ*

^{T}/

*d*

*t*about the tail axis,

*e*^{T}, which is fixed in relation to the cell body, such that

*e*^{T}=−

*e*^{H}. Writing the relative position of a point from the junction as , the surface velocity distribution can be expressed as 2.4 where

*H*denotes the surface of the cell body and

*T*denotes the surface of the flagellum. Some additional equations are needed to determine the unknown velocities,

**,**

*U*

*Ω*^{H}and

*Ω*

^{T}.

Assuming that the bacterium is swimming in the absence of external forces, such as gravity (valid for neutrally buoyant cells), the total force acting on the combined surfaces of the cell body and flagellum is zero. Similarly, no net torque acts on the bacterium. However, there will be a net torque acting on the flagellum balanced out by an opposite torque acting on the body. Two of the torque components cannot be predetermined, but the torque in the direction of the tail axis can be specified (assuming this is mediated by the motor), leading to an additional torque balance equation. These conditions are expressed mathematically as
2.5and
2.6
where *e*^{T} is the unit direction vector of the tail axis and *τ*^{M} is the magnitude of the applied motor torque. In most simulations, we will prescribe a constant value to *τ*^{M}, as bacterial motors produce constant torque over a large range of frequencies. It is also possible to fix the rotation rate, *Ω*^{T}, rather than the torque, as was commonly done in past investigations.

After discretization, we have 3*N*_{N} unknown force components on the mesh nodes, where *N*_{N} is the number of nodes on our surface mesh, and a further seven unknown velocities (** U**,

*Ω*^{H}and

*Ω*

^{T}). The boundary-integral equation (2.2) leads to 3

*N*

_{N}linear equations in these unknowns. Force and torque balance considerations (2.5) give six further equations and the prescription of motor torque (2.6) provides one more. We obtain a dense (3

*N*

_{N}+7)-dimensional linear algebraic problem, which is solved using standard LU decomposition subroutines (Press

*et al.*1996).

### (d) Tracking bacterial trajectories

We implement a method similar to Smith *et al.* (2009) based on Heun’s second-order predictor–corrector algorithm to perform time steps and track the motion of the bacterium. The predictor step involves calculating the instantaneous velocities, *U*_{n}, and at time step *n*. The tentative positions at time step (*n*+1) are found by translating all nodes by *δ*** x**=

*U**δ*

*t*, rotating about the junction by

*δ*

**=**

*θ*

*Ω*^{H}

*δ*

*t*and then rotating the flagellum through an angle

*δ*

*ϕ*

^{T}=

*Ω*

^{T}

*δ*

*t*about its axis. The velocities at the next time step, , and are then calculated assuming this tentative configuration, and the corrected velocities are defined by , and . The corrected positions are then found as before, using the corrected velocities.

This time-stepping process is necessary to accurately predict the motion of the model swimmer, but long-time-scale simulations require significant computational power. We now describe more efficient alternatives for analysing the approximate dynamics.

### (e) Phase-averaged velocities

The bacterium experiences fluctuations in velocity as the flagellum progresses through each revolution. Rather than considering an instantaneous velocity, which depends on the phase of the flagellum, it is often more useful to find the net motion over a cycle. To do this, we take an average from several computations using the same approach as Ramia *et al.* (1993). For this, the tail rotation rate *Ω*^{T} relative to the head is fixed and the arithmetic mean of six instantaneous velocities is calculated for the same head position but varying flagellum phases. This will be denoted . Other phase-averaged quantities, such as the cell-body rotation and motor torque , may be computed similarly, and are used in our analyses in place of their instantaneous values. If the kinematics for constant motor torque are required, the computed average velocity can simply be scaled by the average torque experienced under constant tail-rotation rate, since the swimming speed depends linearly on *Ω*^{T}, which depends linearly on the applied torque.

### (f) Swimming speed and optimization

Variations in power efficiency with the geometrical parameters dictating the number of turns on the flagellum, *N*_{λ}, helix pitch angle, *a**k*, starting region length *k*/*k*_{E}, flagellum radius , flagellum length and cell-body aspect ratio have been described in the past (Higdon 1979*b*; Phan-Thien *et al.* 1987; Fujita & Kawai 2001; see table 1 for a description of parameters). Phan-Thien *et al.* and Fujita & Kawai started with an initial choice of parameters and successively optimized over each parameter in turn, which is essentially one iteration of a direction set method for optimization. We aim to extend the coverage of parameter space to find an improved estimate for the optimal parameters. We also consider torque efficiency in the same way and compare results.

The power efficiency is defined as a ratio of the power needed to push a sphere of the same volume as the cell to the required mechanical power of the swimming motion giving the same mean speed, i.e.
2.7
where is the average value of the power, *P*=*τ*^{M}*Ω*^{T}, required to produce the motion and is the radius of a sphere of the same volume as the bacterium’s head. As the definition indicates, power efficiency is a measure of the square of the swimming speed that can be attained per unit power consumption.

Motivated by the observation that bacterial motors produce relatively constant torques over a wide range of frequencies, we can also look at torque efficiency, defined as the ratio of the torque required to rotate a sphere to the motor torque *τ*^{M} required to propel the organism with mean speed equal to the equatorial surface speed of the sphere, i.e.
2.8
This gives a measure of the swimming speed attained per unit torque exerted by the motor.

Consider the problem of optimization over *N* geometrical parameters. The strategy is to first evaluate the power or torque efficiency on a regular grid in parameter space using the methods outlined in §2*c* and §2*e*. The Matlab function `interpn` is used to produce an *N*-dimensional cubic interpolant of the sampled efficiency data, which is then maximized with `fmincon` within the bounds of the sampled parameter range.

### (g) Phase-plane analysis

The trajectory of a swimming bacterium will clearly depend on its initial state. It would be impractical to perform long simulations for enough different initial configurations to claim an understanding of the complete picture. It is not even clear *a priori* how long each simulation must be to show the long-term behaviour. Instead, we simplify the analysis by examining just a few quantities of interest. The phase-averaged velocity removes dependence on *ϕ*^{T} so it is only necessary to specify two components of the bacterium’s state: the height *h* of the junction position above the no-slip plane boundary and the angle of inclination *θ* of the bacterium with respect to this plane (see figure 2). By symmetry, translations by a vector parallel to the plane and rotations about the plane’s normal direction only change the perspective of the observer. For a spheroidal cell body and coaxial flagellum, rotations about the axis *e*^{H} also leave the configuration unchanged. Hence, we can express the process of boundary-element modelling and phase-averaging symbolically with functions *f* and *g*, leading to the autonomous system
2.9
Here, (assuming, without loss of generality, that the *x*_{2}-direction is the relevant axis of rotation for inclination, corresponding to *e*_{2} in figure 2) and are the phase-averaged rates of change of inclination and height, respectively. By sampling sufficiently many points in the phase plane, we can predict the long-term behaviour starting from any configuration. More specifically, we can predict whether the bacterium will collide with the wall, swim at some steady height above the wall, enter a periodic limit cycle of dips and climbs, or simply swim away from the wall. There may be regions of initial states that give rise to different long-term outcomes and this would also be evident in the phase-plane analysis. For such analyses, Matlab’s built-in functions are used to produce quiver plots of the phase plane and find path lines. Results are given in §3*c*.

Using this method, a lot of computation time is saved as fewer instances need to be calculated. However, only an approximate agreement should be expected between the predicted trajectory and that computed by tracking, since the phase-averaging removes details of the short-time-scale fluctuations, which inevitably perturb the exact path. Therefore, some tracking simulations are carried out to support and verify results derived from phase-plane analysis.

### (h) Locating stable configurations

Trials with phase-plane analysis sometimes indicated the presence of stable configurations, i.e. (*θ**,*h**) satisfying with locally attractive behaviour. The values for (*θ**,*h**) can be estimated from phase-plane analysis, but a quicker (quicker than collecting sufficient phase-plane data) and more accurate alternative is to compute the fixed point directly with a two-dimensional root-finding algorithm. We use Newton’s method with backtracking (Press *et al.* 1996) for this. Once the stable point is obtained, we can calculate the radius of curvature of the circular orbit from the formula
2.10
where is the rotational velocity parallel to the wall, but in the vertical plane containing the bacterial axis, and is the rotational velocity in the direction normal to the wall (directions *e*_{1} and *e*_{3}, respectively, in figure 2). This formula assumes that motion is purely parallel to the wall and is therefore only valid after the swimmer has entered a stable orbit.

In order to evaluate the effects of geometrical parameters on boundary accumulation, we track the dependence of *h** on each parameter. For small parameter adjustments, small changes in *h** and *θ** are expected, so the current steady-state configuration is a good initial estimate for the next run. Typically, two or three iterations of the Newton method are sufficient to accurately locate successive stable configurations as a parameter is varied in small steps.

This scheme breaks down when the parameter passes some critical values. This is either because the steady configuration approaches the boundary (causing the bacterium to descend into the wall in tracking simulations) or the stable height *h** increases asympotically to infinity (representing escape from the surface), and very small steps are required to follow the rapid changes.

A simple method for exploring parameter space is as follows. A coarse grid is constructed to cover the parameter region of interest, which should contain a point for which a stable configuration can be estimated from known data. Starting at this estimate, the stable configuration is obtained and then the stable configurations at each neighbouring grid node are found by branching out. If the destination node is successfully reached and the corresponding stable configuration is found, this is a ‘stable grid node’. For each stable grid node found, the algorithm branches out towards all neighbouring nodes that have not already been found to be stable. If this process fails at any stage, the grid is refined near the point of failure to capture more detail about the boundary of the stable region.

A threshold of (35 nm in dimensional units) is imposed for the separation between the wall and bacterium. Below this, the bacterium is considered to have descended into a region characterized by wall–cell surface nanoscale forces and the possibility of adhesion; this is not included in our model. The length scale characterizing when these surface interactions are important is highly variable, depending on the solutes and their concentrations as well as the cell and wall surface properties. For experiments mimicking physiological conditions, wall–cell surface nanoscale forces are observed to become important for bacteria at separations of 10–60 nm (Klein *et al.* 2003). Thus, in the presentation of our results, we have considered the average, 35 nm, as the cut-off after which we no longer track the bacterium’s descent into the wall. This distance is also consistent with estimates by Vigeant *et al.* (2002).

This procedure results in a representative set of parameter values for which a finite stable separation height exists, and it would be possible to construct the outline that divides the stable from the unstable region in parameter space. Results are given for a two-parameter implementation of this method.

## 3. Results

### (a) Validation of numerical code

Our code was validated by comparing simulations of a sphere translating in free space and falling towards a no-slip boundary with exact solutions (Brenner 1961). We also tested a wide range of discretization levels for both the cell body and the flagellum to ensure that results were accurate to about 1 per cent compared with much finer discretizations. Further details of testing are given in §3 of the electronic supplementary material.

Simulation results were in excellent agreement with those published by other authors such as Phan-Thien *et al.* (1987) and Higdon (1979*b*). Our results tended to be slightly closer to Higdon’s than to Phan-Thien’s. A comparison is shown in figure 3.

### (b) Optimization

Power efficiency appears to increase monotonically with decreasing flagellum radius (Phan-Thien *et al.* 1987), so we will fix this parameter at a value , where is the length scale used for non-dimensionalization (refer to table 1 for a description of parameters). This is the right order of magnitude, though larger than observed in most monotrichous bacteria (where –0.04, estimated from data in Brennen & Winet 1977), which allows greater computational accuracy without excessive mesh refinement. In addition, we will fix the value *k*/*k*_{E}=1 since efficiency is insensitive to the starting region length within a realistic range applicable to a bacterial model. Finally, only axisymmetric spheroids will be considered, i.e. , since this is a good approximation for rod-shaped bacteria. Thus, we will optimize swimming efficiency over four parameters: number of turns on flagellum *N*_{λ}, pitch angle *a**k*, flagellar length and cell-body aspect ratio .

To begin with, consider the problem of optimization in free space over two parameters: wavelength, , and flagellar length, . This is equivalent to optimizing over *N*_{λ} and since these parameters are inter-related. Figure 4*a* shows the variation in power efficiency with these two parameters when we fix *a**k*=1 and . We find that efficiency is within about 10 per cent of the maximum value in the region , .

The results in figure 4*b* are for bacteria oriented parallel to the boundary at a distance , which is found to be a typical stable swimming height (see §3*d*). The overall efficiency is reduced slightly near the boundary, but the location of the optimum does not change significantly. Even upon inclusion of all four geometrical parameters, the configurations optimized for free-space swimming and those optimized for boundary swimming do not differ greatly and are almost equally efficient.

However, the optimum shape is quite different if we consider torque efficiency, *η*^{τ}. As figure 4*c* shows, there is a much greater preference for short, small-wavelength helices. Considering only two parameters, the optimal region is roughly and . Again, the presence of the boundary does not significantly affect the optimum in either the two- or four-parameter case. Results of geometrical optimization in each case are listed in table 2 and comparisons with experimental measurements are given in table 3. The power-optimal parameters are very close to those previously reported for free-space swimming (Phan-Thien *et al.* 1987).

### (c) Phase-plane analysis of cell dynamics above a no-slip boundary

As one would expect, bacteria tilted too heavily towards the wall will descend into the influence of the wall–cell surface nanoscale interactions, while those initially tilted too much away from the wall will escape. Apart from such extremes, the final state is insensitive to initial conditions and depends only on the bacterial geometry. Some bacteria are deflected and swim away from the wall, others turn and swim towards the wall while the rest settle at a stable configuration. As described in previous analyses (Lauga *et al.* 2006), the swimmer must rotate about the normal direction to the boundary in order to balance the torque induced by the body and tail rotating axially near a boundary. This leads to a circular orbit above the wall.

Figure 5 is an example of a phase-plane diagram that we find by evaluating phase-averaged velocities of a bacterium as described in §2*g*. Here, there is a stable fixed point in the phase plane corresponding to the bacterium swimming (on average) at a constant height above the wall and constant tilt. This result was verified by performing a tracking simulation (also shown in figure 5).

Generally, the predicted stable points are very close to the stable configuration manifested in tracking simulations. For example, in figure 5, the predicted stable point is at , while the simulated trajectory oscillates around .

### (d) Boundary accumulation and variation of parameters

While stable spirals are often found, some bacterial geometries have phase portraits in which all trajectories either approach or escape from the wall. To understand why certain bacteria have a tendency to escape while others maintain a stable separation from the wall, the variations in stable height with some key geometrical parameters are plotted in figure 6. The stable inclination angle *θ** was found to be positive (body pointing away from wall) in all tested cases and was never more than about 5^{°} from horizontal. Since inclination is not an indicator of boundary accumulation, it will not be discussed in detail here.

We see from figure 6 that longer flagella, higher numbers of turns and thicker filaments result in closer accumulation heights, while tightly wound or high-amplitude helices (large *a**k*) and long starting region lengths push the swimmer further away from the wall. Body shape is also important, with elongated cells swimming further from the wall than more oblate cells. Apart from variations in helical pitch angle, the radius of curvature of the stable circular orbit tends to follow the same trend as the accumulation height, i.e. the orbit radius generally increases with . However, figure 6*c* shows that the radius achieves a maximum near *a**k*=1, despite the fact that the stable height, *h**, increases monotonically with *a**k* throughout the tested range.

Figure 7 indicates a region of parameter space where boundary accumulation is found to occur when the cell-body aspect ratio and the length of the flagellum are allowed to vary. For most stable configurations, it was found that the accumulation height, *h**, was in the range , where is the radius of a sphere with the volume of the cell body. The wall separation distance is roughly less than the stable height, which is measured from the cell body–flagellum junction. For a fixed flagellum length, increasing the aspect ratio of the cell body, , results in an increased accumulation height. Above a critical value, which varies with model parameters, the increase in *h** becomes very rapid and continues to the point where the cell may be considered to have escaped wall effects. Decreasing the aspect ratio lowers the accumulation height until the bacterium descends into the boundary.

The accumulation height generally increases as the flagellum length is reduced, although there can be fluctuations as well. This is most noticeable at low aspect ratios, leading to the complicated boundary contour on the left-hand side of figure 7, but is also visible in figure 6*d*. In figure 7, we find that, at aspect ratios greater than 1.4, the accumulation height increases very rapidly below a critical flagellar length, which varies with the aspect ratio, and the cell effectively escapes from the wall. Figure 8 shows the same regions of surface accumulation behaviour laid over plots of power and torque efficiency. Figure 9 plots trajectories illustrating the three types of long-term behaviour—descent, stable orbiting and escape from the wall—as predicted in figure 7.

## 4. Discussion and conclusions

### (a) Optimization

The plots in figure 4 show rather broad optima in power and torque efficiencies when considering variations in flagellum length and wavelength. Configurations achieving at least 90 per cent of the optimal efficiency cover parameter ranges spanning a factor of two. This suggests that biological fine tuning is not required for efficient swimming; hence, efficient swimming is robust to biological variation. There is also some freedom for the cell to deviate from the swimming optimum in order to accommodate other criteria.

Both power and torque efficiency are relatively insensitive to the presence of the boundary, at least up to separations of the order of . Thus, it is likely that no special adaptation is required for efficient swimming near surfaces. In other words, efficient free-space swimmers are also efficient near boundaries. This is consistent with the observation that parallel swimming speeds close to a surface are indistinguishable from swimming speeds far from a surface (Frymier & Ford 1997).

Given the marked difference between power-optimal and torque-optimal shapes, we briefly consider the extent to which observed bacterial shapes and speeds resemble the predicted power and torque optima. Flagellar bundle frequencies for *Escherichia coli* are estimated to be around 100–200 Hz, depending on experimental conditions (Darnton *et al.* 2007). For a frequency of 130 Hz, the monotrichous power-optimized configuration is predicted to swim at a speed of about 10 body lengths per second. This agrees roughly with observations for bacterial swimming in media with viscosity similar to water (Darnton *et al.* 2007), whereas torque-optimized shapes are predicted to require rotation rates 10–15 times higher to achieve the same swimming speed. However, in *E. coli* and a number of other species, the constant torque regime of motors holds only up to speeds of around 100–200 Hz (Sowa & Berry 2008). Hence, we observe that constant torque optima require rotation rates that are beyond the approximately constant torque region of the torque–speed curve for bacterial motors (except in sufficiently viscous media) if we are to achieve commonly measured swimming speeds.

Table 3 illustrates the approximate shapes of several bacterial species, as well as our predicted power and torque optima. Note that there is considerable variation between species; *Photobacterium phosphoreum* closely resembles the power optimum and *V. cholera* is reminiscent of the torque optimum, but other bacteria have characteristics between or unlike either choice of optima. For example, some bacteria have as many as four turns on their flagella, while all of our computations for optimal shapes predict 1.5 or fewer turns. Thus, we observe that elements associated with power and torque efficiency can be seen in some monotrichous bacteria, but there is certainly no indication of universal efficiency principles.

More generally, the results show that different, but plausible, objective functions can give very different results in optimization. It is therefore necessary to justify the choice of objective, as well as possible constraints in any optimization study.

### (b) Boundary accumulation

Consistent with the simulation study of spermatozoa by Smith *et al.* (2009), we find that some, but not all, flagellated microswimmer shapes can lead to accumulation at a finite separation from a wall, purely as a result of hydrodynamic forces. When boundary accumulation is present, we observe a small, positive stable tilt angle, *θ**, also consistent with the previous study. Although a detailed investigation has not been carried out, we suggest that pointing slightly away from the wall is required to balance the attractive effect of near-wall swimming discussed by Berke *et al.* (2008). Smith & Blake (2009) describe the same phenomenon in more detail for simulated spermatozoa.

Accumulation separations predicted for spermatozoa (Smith *et al.* 2009) were approximately half a flagellum length (), while the present study finds typical separations of about half a head radius () for bacteria. This is in line with the estimate of 0–2 μm for *E. coli* by Frymier *et al.* (1995), but about an order of magnitude larger than more recent measurements by Vigeant *et al.* (2002). It is worth noting that both of these experimental studies used *E. coli*, and it is not clear how the dynamics might differ for a flagellar bundle compared with a single flagellum.

The phase-plane plots suggest that the stable configuration has a large basin of attraction, whereas the previous work on spermatozoa (Smith *et al.* 2009) reported a somewhat sensitive dependence on the initial tilt angle, and observed that the same cell could swim steadily parallel to the surface or ‘bounce’ off and escape. This phenomenon has not been seen in the current bacterial model, which may be accounted for by differences between spermatozoa and bacteria. Not only are flagella of spermatozoa much longer and thicker in proportion to head size, but the beat patterns are also different.

It was noted that previous simulations of bacteria swimming near solid walls have generally descended into the wall rather than settling at some stable configuration. This can be explained in light of the current findings, such as those shown in figure 7. Cell bodies with low aspect ratio, such as the commonly used sphere, have a strong tendency to descend into the boundary. It would be interesting to ascertain experimentally whether spherical bacteria indeed collide more often or remain closer to surfaces than rod-shaped specimens.

This finding has implications for the representation of microswimmers by a small collection of singularities. There have been studies that treat the organism as a Stokes dipole, representing a thrust from the flagellum and a drag from the cell body (e.g. Wolgemuth 2008). This kind of simplification gives reasonable far-field behaviour and allows scaling up to investigate the dynamics of large numbers of swimmers. However, it is not obvious how to form a simple representation that still captures the near-field complexities of flagellar motility (Smith & Blake 2009). The results described here reinforce the fact that bacterial shape significantly affects locomotion near boundaries and any accurate representation would have to accommodate these effects. Other numerical techniques, such as the method of regularized Stokeslets, have been used to analyse flow patterns of multiple swimmers (Cisneros *et al.* 2007), and it would be interesting to compare the results of these different methods.

figure 6*a*,*c* suggests that conformational changes in the flagellum could affect boundary-accumulation behaviour. It is known that the flagella of *R. sphaeroides* can take the normal helical form during swimming or a high-amplitude, short-wavelength form typically when the motor stops or runs at low frequency (Armitage *et al.* 1999). If multiple conformations can be used for propulsion, then results indicate the possibility that some conformations lead to accumulation, while others lead to escape from surfaces.

The radius of curvature of circular trajectories near surfaces is generally between 10 and 50 μm (summarized by Lauga *et al.* 2006). Using a physiological range of parameter values, we predict radii of curvature between 10 and 300 μm, heavily dependent on certain parameters such as the aspect ratio of the cell body. Different flagellar conformations would therefore lead to different path curvatures. If flagellar conformations during forward and backward swimming are sufficiently different, this effect may contribute to the factor of 10–20 difference in path curvatures found by Kudo *et al.* (2005) between these two modes of swimming in *V. alginolyticus* cells near surfaces.

The influence of cell size and aspect ratio on the trajectories of bacteria near surfaces has been observed before. Hulme *et al.* (2008) used the correlation between cell length and trajectory curvature to design a device that could sort populations of *E. coli* by cell length. As the cells mature, they become longer without growing wider. In other words, the aspect ratio increases, and the present model predicts that this will increase the separation distance from the channel walls and the radius of curvature of the resulting circular path, which is consistent with observation. By applying techniques described in this work to specific bacterial species, it would be possible to improve the design of sorting devices. However, it should be noted that the present model has only been developed for a half-space domain and would require modifications to accurately describe motion in narrow channels.

One of the key questions to ask at this point is how the requirements for boundary accumulation relate to power and torque efficiency. As shown in figure 8, the region of high-power efficiency coincides roughly with part of the region of stable accumulation behaviour when considering variations in the cell-body aspect ratio and flagellum length. Within the investigated range of these parameters, all boundary-accumulating geometries were fairly power efficient (greater than 70% of optimal power efficiency). In particular, the power optimum exhibits boundary accumulation. However, for all but the highest power efficiencies, it was possible to find geometries that descend into the wall and others that escape from the wall. In other words, power-efficient cells can be formed with any of the three near-surface swimming behaviours. In contrast, high torque efficiency requires shorter flagella, which result in escape from the boundary. Thus boundary-accumulating swimmers do not achieve optimal torque efficiency, though some can reach up to 85 per cent of the optimum.

## 5. Summary

We have developed a boundary-element model for a flagellated bacterium swimming in a quiescent Newtonian fluid and investigated the motion resulting from rotation of the helical flagellum when the fluid is unbounded, as well as when the bacterium is near a no-slip plane boundary. We have found body and flagellum shapes that give optimal swimming efficiency in terms of power and torque. In the case of power, this optimum agrees well with estimates from previous simulation studies. Torque efficiency has not been considered before, yet it is biologically interesting given that bacterial motors produce nearly constant torque over a wide range of operating frequencies. High torque efficiency means fast swimming under this constraint. Torque and power efficiency considerations lead to different optima; hence optimization studies clearly need to justify optimization objectives.

It is probable that a combination of factors determines biological fitness such that most flagellated bacteria are of an intermediate design. Nonetheless, the presence of broad optima suggests that biological fine tuning is not crucial for efficient swimming. In addition, it was determined that, despite the dynamical influence of no-slip surfaces, proximity to such a surface does not significantly alter the optimal shape for swimming whether power or torque efficiency is considered.

This is the first in-depth numerical study of how bacteria morphology affects boundary accumulation. We have shown that details of the bacterium’s shape have a profound influence on the tendency for the cells to accumulate close to solid surfaces; bacteria can be hydrodynamically attracted or repelled by walls simply by modifying geometrical parameters of the body and flagellum. Using phase-plane analysis techniques that have not been applied to microswimmer dynamics before, we predict that bacteria with long, thick flagella and short bodies are drawn closer to walls than those with short, thin flagella and elongated bodies. For most of the parameter space explored, accumulation heights were between and , where is the radius of a sphere with the volume of the cell body. Furthermore, cells settle at a configuration with the body tilted a few degrees away from the wall.

## Acknowledgements

This publication is based on work supported in part by award no. KUK-C1-013-04, made by King Abdullah University of Science and Technology (KAUST). D.J.S. acknowledges the MRC (Special Training Fellowship G0600178 in Computational Biology). The authors would like to thank Prof. John Blake of the School of Mathematics, University of Birmingham, Dr Jackson Kirkman-Brown of the Centre for Human Reproductive Science, Birmingham Women’s Hospital and Mr Hermes Gadêlha of the Mathematical Institute, University of Oxford for invaluable discussions.

## Footnotes

- Received October 4, 2009.
- Accepted December 1, 2009.

- © 2010 The Royal Society