Using Experimentally Calibrated Regularized Stokeslets to Assess Bacterial Flagellar Motility Near a Surface

Abstract

The presence of a nearby boundary is likely to be important in the life cycle and evolution of motile flagellate bacteria. This has led many authors to employ numerical simulations to model near-surface bacterial motion and compute hydrodynamic boundary effects. A common choice has been the method of images for regularized Stokeslets (MIRS); however, the method requires discretization sizes and regularization parameters that are not specified by any theory. To determine appropriate regularization parameters for given discretization choices in MIRS, we conducted dynamically similar macroscopic experiments and fit the simulations to the data. In the experiments, we measured the torque on cylinders and helices of different wavelengths as they rotated in a viscous fluid at various distances to a boundary. We found that differences between experiments and optimized simulations were less than 5% when using surface discretizations for cylinders and centerline discretizations for helices. Having determined optimal regularization parameters, we used MIRS to simulate an idealized free-swimming bacterium constructed of a cylindrical cell body and a helical flagellum moving near a boundary. We assessed the swimming performance of many bacterial morphologies by computing swimming speed, motor rotation rate, Purcell’s propulsive efficiency, energy cost per swimming distance, and a new metabolic energy cost defined to be the energy cost per body mass per swimming distance. All five measures predicted that the optimal flagellar wavelength is eight times the helical radius independently of body size and surface proximity. Although the measures disagreed on the optimal body size, they all predicted that body size is an important factor in the energy cost of bacterial motility near and far from a surface.

1. Introduction

Living organisms emerge, evolve, and reside within habitats, and the physical interactions among organisms and their environments impose selective forces on their evolution. In their low Reynolds number surroundings, bacteria such as Escherichia coli and Pseudomonas aeruginosa have evolved a mechanical motility system to propel themselves through fluids. This system consists of one or more helical flagella, and these flagellar organelles are attached to the cell body by rotary nanomotors. Flagellar motor rotation is driven by an ion flow through the motor, causing the flagellum and the bacterial cell body to rotate in opposite directions [1]. A bacterium swimming through a fluid can be described as a non-inertial system in which the mechanical power output by the motor is instantaneously dissipated by fluid drag on the body and flagellar filaments. The interaction between the bacterium and the fluid generates a flow that results in the net motion of the bacterium. Different flows can be more or less favorable to the survival of an organism [2], and the presence of a surface introduces boundary effects that modify how a swimming cell interacts with the fluid [3]. We consider here the example of a unicellular motile flagellate bacterium swimming through a fluid near to a surface and how the conformation of the bacterial cell body and the flagellar organelle may be optimized for such an environment; see Figure 1 for an image of our model.

The efficiency of the bacterial motility system has been the focus of numerous theoretical [4,5,6], computational [7,8,9,10,11,12], and experimental works [13,14,15]. For a comprehensive review, see Ref. [2]. In an early paper on swimming efficiency, E. Purcell discussed two measures: the propulsive efficiency (Purcell efficiency) and the energy consumed during bacterial motion per body mass [4]. The Purcell efficiency—a specialized form of the Lighthill efficiency [16] for rotary motor-driven bacterial propulsion—is defined as the ratio of the least power needed to translate a bacterial body against fluid drag to the total power output by the motor during motion of the bacterium. Most work has focused on the Purcell efficiency because it is a scale-independent function of the geometries of the cell body and flagellum. One shortcoming of this measure, however, is that it is independent of the motor’s response to an external load imposed by the environment and therefore cannot assess the biological fitness of the bacterial motor. Another measure of bacterial performance used by a few authors is the distance traveled by a bacterium per energy input by the motor [14,17], which provides a different means of evaluating fitness, as explained below.

In this work, we investigate and compare predictions of the optimal bacterial motility system made by five measures. The first two measures are related directly to the motion of a bacterium: the swimming speed and the motor rotation frequency. Swimming speed is important because bacteria live in an environment where nutrients diffuse on time and length scales comparable to bacterial motion. To effectively achieve chemotaxis, bacteria must move quickly enough to sample their chemical environment before it is randomized by diffusion [4,12]. Motor rotation rate is important because the bacterial motor has a characteristic frequency response that depends on the external torque load [18,19,20,21].

The other three performance measures we studied are based on the mechanical energy cost to achieve motility: the Purcell inefficiency (or the inverse of the Purcell efficiency), the inverse of distance traveled per energy input, and the metabolic energy cost, which we define to be the energy output by the motor per body mass per distance traveled. Each of these measures compares the ratio of the power output of the bacterial motor to the performance of a particular task. The rationale for introducing the metabolic cost function is that it measures the actual energetic cost to the organism to perform a specific biologically relevant task, i.e., translation through the fluid. Additionally, both the energy consumed per distance traveled and the metabolic energy cost depend upon the rotation speed of the motor. Thus, their predictions about optimal morphologies depend upon the torque–speed response of the motor.

To determine the values of performance measures attained by different bacterial geometries, we employed the method of regularized Stokeslets (MRS) [22] and the method of images for regularized Stokeslets (MIRS) [23], the latter of which includes the effect of a solid boundary. Employing MRS and MIRS requires determining values for two kinds of free parameters: those associated with computation and those associated with the biological system. As with any computational method, the bacterial structure in the simulation is represented as a set of discrete points. The body forces acting at those points are expressed as a vector force multiplied by a regularized distribution function, whose width is specified by a regularization parameter. Though other simulations have produced numerical values for dynamical quantities such as torque [24] that are within a reasonable range for bacteria, precise numbers are not possible without an accurately calibrated method.

In this work, we present for the first time in the literature a method for calibrating the MIRS using dynamically similar experiments. There is no theory that predicts the relationship between the discretization and regularization parameters, though one benchmarking study showed that MRS simulations could be made to match the results of other numerical methods [25]. To determine the optimal regularization parameter for chosen discretization sizes, we performed dynamically similar macroscopic experiments using the two objects composing our model bacterium: a cylinder and a helix, see Figure 1. Such an approach was previously used to evaluate the accuracy of various computational and theoretical methods for a helix [26], but the study did not consider the effects of a nearby boundary.

By measuring values of the fluid torque acting on rotating cylinders near a boundary, we verified the theory of Jeffery and Onishi [27], which is also a novelty in our work. We then used the theory to calibrate the ratio of discretization to regularization size in MRS and MIRS simulations of rotating cylindrical cell bodies. Because there are no exact analytical results for helices, we determined regularization parameters for helices that were discretized along their centerlines by fitting simulation results directly to experimental measurements. Calibrating our simulations of rotating cylinders and helices with the experiments allowed us to build a bacterial model with a cylindrical cell body and a helical flagellum whose discretization and regularization parameter are optimized for each part.

To impose motion on the bacterial model, we needed only to specify the motor rotation—a consequence of there being no body forces acting on the bacterium [24]. The motor rotation rate, however, depends upon the external load [14,18,19,20]. A novel aspect of our simulation method was to ensure that the motor rotation rate and the torque load matched points on the experimentally determined torque–speed curve [18,21]. The dynamical quantities output from the simulations were then used to compute swimming performance measures for different bacterial geometries at various distances from the boundary. Among these measures, we defined a new metabolic energy cost that quantifies the energy per body mass required for bacterial propulsion, which provides a new tool for analyzing the efficiency of bacterial swimming.

Our paper is organized as follows: Section 2 discusses our implementation of the MRS and the MIRS, our use of dynamically similar experiments to calibrate the simulations, and our determination of the torque–speed response curve for the motor; Section 3 compares our five fitness measures: free swimming speed, motor frequency, inverse Purcell efficiency, energy per distance, and metabolic energy cost; and Section 4 discusses the predictions made by each fitness measure and comments on future directions of our work.

2. Materials and Methods

2.1. Numerical Methods

Bacterial motility using a helical flagellum often involves multiple flagella, and bodies may be spherical, cylindrical, or helical [28]. We reduced the complexity by considering a simpler biomechanical system of a regular cylindrical body to which a single, uniform flagellum is attached, as shown in Figure 1. This simple system, however, contains the same essential geometric factors as bacteria such as E. coli, which have a long rod-shaped body and helical flagella that bundle together, forming a single helix. Our goal was to assess how the performance of our model organism changes when its geometrical parameters and distance to an infinite plane wall are varied in numerical simulations. We quantified the performance of different models by computing speed, motor rotation rate, and the three energy cost measures. A glossary of symbols used in the bacterial models the and the calculated energy measures is displayed as

Table 1. Glossary of parameters for the computational and experimental work.

Dynamic Viscosity of the Fluid	μ	Distance of Flagellum to Wall	d
Cylindrical cell body		Helical flagellum
Geometrical parameters		Geometrical parameters
Length	ℓ	Axial length	L
Radius	r	Helix radius	R
		Wavelength	λ
		Filament radius	a
Computational parameters		Computational parameters
Optimal discretization factor	γc	Optimal filament factor	γf
Regularization parameter	ϵc	Regularization parameter	ϵf
Discretization size	dsc	Discretization size	dsf
Body mass	m	Motor angular frequency	Ωm
Axial drag force	F	Axial torque	τ
Swimming speed	U	Purcell inefficiency	${\displaystyle {\mathcal{E}}_{\mathrm{Purcell}}^{-1}\equiv \frac{\tau {\mathrm{\Omega }}_m}{FU}}$
Energy per distance traveled	${\displaystyle E^{*}\equiv \frac{\tau {\mathrm{\Omega }}_m}{U}}$	Metabolic energy cost	${\displaystyle \mathcal{E}\equiv \frac{E^{*}}{m}}$

Lengths (ℓ, r, L, $\lambda$, a, and d) are made scale-free by dividing by the helical radius R. See Figure 1 for image of the model.

.

We composed our model of a bacterium with a cylindrical cell body and a tapered left-handed helical flagellum as shown in Figure 1 and Figure 2. The flagellar centerline is described by

$$\left\{\begin{array}{cc}\hfill x\left(s\right)& =(1-e^{-k^2{s}^2})Rsin(ks+\theta )\hfill \\ \hfill y\left(s\right)& =(1-e^{-k^2{s}^2})Rcos(ks+\theta )\hfill \\ \hfill z\left(s\right)& =s\hfill \end{array}\right.$$

(1)

where $0\le s\le L$ with L the axial length in the z-direction, k is the wavenumber $2\pi /\lambda$ with $\lambda$ the wavelength, and $\theta$ is the phase angle of the helical flagellum at 16 evenly spaced phases.

The parameter values used for the bacterium models shown in Figure 2 are given in

Table 2. Parameters values used in numerical simulations.

Parameter	Value	Unit	Reference
μ	0.93	10−3 Pa·s
Cell body
ℓ	(a)	$\mathsf{\mu }\mathrm{m}$	[21]
r	(b)	$\mathsf{\mu }\mathrm{m}$	[21]
${\gamma }_c$	$6.4$
$d{s}_c$	$0.096$	$\mathsf{\mu }\mathrm{m}$
${ϵ}_c=d{s}_c/{\gamma }_c$	$0.015$	$\mathsf{\mu }\mathrm{m}$
Flagellum
L	$8.3$	$\mathsf{\mu }\mathrm{m}$	[21]
$\lambda$	(c)	$\mathsf{\mu }\mathrm{m}$
R	$0.2$	$\mathsf{\mu }\mathrm{m}$	[21]
a	$0.012$	$\mathsf{\mu }\mathrm{m}$	[21]
${\mathrm{\Omega }}_m/\left(2\pi \right)$	154	$\mathrm{Hz}$	[21]
${\gamma }_f$	$2.139$
$d{s}_f={ϵ}_f$	$0.026$	$\mathsf{\mu }\mathrm{m}$
d	(d)	$\mathsf{\mu }\mathrm{m}$

(a) $\ell \in$ {1.9, 2.2, 2.5, 2.8, 3.1} $\left(\mathsf{\mu }\mathrm{m}\right)$. (b) $r\in$ {0.395, 0.4175, 0.44, 0.4625, 0.485} $\left(\mathsf{\mu }\mathrm{m}\right)$. (c) $\lambda \in$ {0.2, 0.5, 0.8, 1.1, 1.4, 1.7, 2.02, 2.22, 2.3, 2.42, 2.6, 2.9, 3.2, 3.6, 4.0, 5.0, 7.0, 9.0 } $\left(\mathsf{\mu }\mathrm{m}\right)$. (d) $d\in$ {0.55, 0.62, 0.71, 0.82, 0.96, 1.12, 1.32, 1.56, 1.85, 2.20, 2.26, 2.52, 2.81, 3.14, 3.5, 3.93, 4.4, 4.93, 5.53, 6.2, 8.2, 10.2} $\left(\mathsf{\mu }\mathrm{m}\right)$.

.

2.1.1. Method of Regularized Stokeslets

The microscopic length and velocity scales of bacteria ensure that fluid motion at that scale can be described using the incompressible Stokes equations. We used the MRS in three dimensions [22] to compute the fluid–bacterium interactions due to the rotating flagellum in free space at steady state:

$$\begin{array}{cc}\hfill \mu ▵\mathbf{u}\left(\mathbf{x}\right)-\nabla p\left(\mathbf{x}\right)& =-\mathbf{F}\left(\mathbf{x}\right)\hfill \\ \hfill \nabla ·\mathbf{u}\left(\mathbf{x}\right)& =0\hfill \end{array}$$

(2)

$\mathbf{u}$ is the fluid velocity, p is the fluid pressure, and $\mu$ is the dynamic viscosity. $\mathbf{F}$ is the body force represented as ${\mathbf{f}}_k{\varphi }_{ϵ}(\mathbf{x}-{\mathbf{x}}_k)$, where ${\mathbf{f}}_k$ is a point force at a discretized point ${\mathbf{x}}_k$ of the bacterium model. In our simulations, we used the blob function ${\displaystyle {\varphi }_{ϵ}(\mathbf{x}-{\mathbf{x}}_k)=\frac{15{ϵ}^4}{8\pi {(r_k^2+{ϵ}^2)}^{\frac{7}{2}}}}$, where $r_k=∥\mathbf{x}-{\mathbf{x}}_k∥$. This radially symmetric smooth function depends on a regularization parameter $ϵ$ which controls the spread of the point force ${\mathbf{f}}_k$. Given N such forces, the resulting velocity at any point $\mathbf{x}$ in the fluid can be computed as

$$\mathbf{u}\left(\mathbf{x}\right)=\frac{1}{8\pi \mu }\sum _{k=1}^N\frac{{\mathbf{f}}_k(r_k^2+2{ϵ}^2)}{{(r_k^2+{ϵ}^2)}^{\frac{3}{2}}}+\frac{({\mathbf{f}}_k·(\mathbf{x}-{\mathbf{x}}_k))(\mathbf{x}-{\mathbf{x}}_k)}{{(r_k^2+{ϵ}^2)}^{\frac{3}{2}}}=\frac{1}{8\pi \mu }\sum _{k=1}^N{S}_{ϵ}(\mathbf{x},{\mathbf{x}}_k){\mathbf{f}}_k$$

(3)

Evaluating Equation (3) N times, once for each ${\mathbf{x}}_k$, yields a $3N\times 3N$ linear system of equations for the velocities of the model points. In the limit as $ϵ$ approaches 0, the resulting velocity $\mathbf{u}$ approaches the classical singular Stokeslet solution. In practice, the specific choice of $ϵ$ may depend on the discretization or the physical thickness of the structure.

In our bacterium model, we discretized the cell body as $N_c$ points on the surface of a cylinder, and we modeled the flagellum as $N_f$ points distributed uniformly along the arc length of the centerline. In Section 3, we present the optimal regularization parameter for the cylindrical cell we obtained by calibrating the simulations based on the experiments and theory. The regularization parameter for the helical flagellum was found by calibrating simulations with experiments, since there is no exact theory for rotating helices, as presented in Section 2.3.

2.1.2. Method of Images for Regularized Stokeslets

We used the method of images for regularized Stokeslets [23] to solve the incompressible Stokes equations (Equation (2)) and simulate bacterial motility near a surface. In the method, the no-slip boundary condition on an infinite plane wall is satisfied by imposing a combination of a Stokeslet, a Stokeslet doublet, a potential dipole, and rotlets at the image point ${\mathbf{x}}_k^{*}$ of each discretized point ${\mathbf{x}}_k$. The image point ${\mathbf{x}}_k^{*}$ is the point obtained by reflecting ${\mathbf{x}}_k$ across the planar surface. The resulting velocity at any point $\mathbf{x}$ in the fluid bounded by a plane can be found in Ref. [23] and is written in the compact form similar to Equation (3):

$$\mathbf{u}\left(\mathbf{x}\right)=\frac{1}{8\pi \mu }\sum _{k=1}^N{S}_{ϵ}^{*}(\mathbf{x},{\mathbf{x}}_k){\mathbf{f}}_k$$

(4)

2.1.3. Force-Free and Torque-Free Models

For a free-swimming bacterium, the only external forces acting are due to the fluid–structure interaction. A bacterium is a non-inertial system, so the net external force and net external torque acting on it must vanish. This means that ${\mathbf{F}}_c+{\mathbf{F}}_f=\mathbf{0}$ and ${\tau }_c+{\tau }_f=\mathbf{0}$, where ${\mathbf{F}}_c$ / ${\tau }_c$ and ${\mathbf{F}}_f$ / ${\tau }_f$ represent, respectively, the net fluid forces and torques acting on the cell body and flagellum. These force-free and torque-free constraints require the cell body and flagellum to counter-rotate relative to each other. In our simulations, the point connecting the cell body and the flagellum ${\mathbf{x}}_r$ represented the motor location, and was used as the reference point for computing torque and angular velocity.

Given an angular velocity ${\mathrm{\Omega }}_m$ of the motor, the relationship between the lab frame angular velocities of the flagellum and the cell body is ${\mathrm{\Omega }}_f={\mathrm{\Omega }}_c+{\mathrm{\Omega }}_m$ [24]. Since ${\mathrm{\Omega }}_m$ is the relative rotational velocity of the flagellum with respect to the cell body, the resulting velocity $\tilde{\mathbf{u}}\left({\mathbf{x}}_k\right)$ at a discretized point ${\mathbf{x}}_k$ on the flagellum ($k=1,\dots ,N_f$) can be computed as ${\mathrm{\Omega }}_m\times {\mathbf{x}}_k$ (this velocity is set to zero at a discretized point on the cell body). Using the MRS (or MIRS) and the six added constraints from the force-free and torque-free conditions, we formed a $(3N+6)\times (3N+6)$ linear system of equations to solve for the translational velocity $\mathbf{U}$ and angular velocity ${\mathrm{\Omega }}_c$ of the cell body and the internal force ${\mathbf{f}}_k$ acting at the discretized point ${\mathbf{x}}_k$ of the model:

$$\begin{array}{cc}\hfill \tilde{\mathbf{u}}\left({\mathbf{x}}_j\right)& =\frac{1}{8\pi \mu }\sum _{k=1}^N{G}_{ϵ}({\mathbf{x}}_j,{\mathbf{x}}_k){\mathbf{f}}_k-\mathbf{U}-{\mathrm{\Omega }}_c\times ({\mathbf{x}}_j-{\mathbf{x}}_r),j=1,\dots ,N\hfill \\ \hfill \sum _{k=1}^N{\mathbf{f}}_k& =\mathbf{0},\sum _{k=1}^N({\mathbf{x}}_k-{\mathbf{x}}_r)\times {\mathbf{f}}_k=\mathbf{0}\hfill \end{array}$$

(5)

where $G_{ϵ}$ is $S_{ϵ}$ from Equation (3) for swimming in a free space or $S_{ϵ}^{*}$ from Equation (4) for swimming near a plane wall. Each ${\mathbf{f}}_k$ represents a point force acting at point ${\mathbf{x}}_k$, which is in principle an internal contact force due to interactions with the points on the bacterium that neighbor ${\mathbf{x}}_k$. Each ${\mathbf{f}}_k$ is balanced by the hydrodynamic drag that arises from a combination of viscous forces and pressure forces exerted on the point ${\mathbf{x}}_k$ by the fluid (Equation (2)).

By computing each ${\mathbf{f}}_k$, we were able to deduce the fluid interaction with each point of the bacterial model. Equation (5) shows that the calculated quantities $\mathbf{U}$, ${\mathrm{\Omega }}_c$, ${\mathbf{F}}_c$, and ${\tau }_c$ depend linearly on the angular velocity ${\mathrm{\Omega }}_m$ since $\tilde{\mathbf{u}}\left({\mathbf{x}}_j\right)={\mathrm{\Omega }}_m\times {\mathbf{x}}_j$.

2.2. Torque–Speed Motor Response Curve

The singly flagellated bacteria we simulated move through their environment by rotating their motor, which causes their body and flagellum to counter-rotate accordingly. Drag force from the fluid exerts equal magnitude torques on the body and the flagellum, and the value of the torque equals the torque load applied to the motor. The relationship between the motor rotation rate and the torque load is characterized by a torque–speed curve, which has been measured experimentally in several organisms [14,18,19,20,21]. In the context of motor response characteristics, speed refers to frequency of rotation. We estimated the torque–speed curve for E. coli with typical values taken from the literature [18,21] to match the body and flagellum parameters also taken from measurements on E. coli [21].

The fluid torque exerted on a rotating object is proportional to its rotation rate under constant environmental conditions in Stokes flow, and thus, plotting the fluid torque versus rotation rate in fixed conditions yields a straight line. Figure 3 shows examples of these ‘load lines’ computed for our bacterial model at different distances from the boundary: the shallower blue line is calculated for a bacterium far from the boundary, and the steeper red line is calculated near the boundary. The load lines shown in Figure 3 were computed with typical body and flagellum parameters for E. coli [21].

The torque–speed curve of the E. coli motor has been determined experimentally by measuring the rotation rate of a bead attached to a flagellar stub and then computing the torque on the bead due to fluid drag. By performing the measurement in fluids of different viscosities, many points on the torque–speed curve were assembled. It was found that the torque–speed curve of the E. coli bacterial motor decreases monotonically from a maximum stall torque (i.e., the zero-speed torque) of about 1300 pN·nm to zero torque, which occurs at a maximum speed of 350 Hz [18,20,21]. There are two linear operating regimes: a low-speed regime from 0–175 Hz and a high-speed regime 175–350 Hz. In the low-speed regime below 175 Hz, the torque is a relatively flat function of the motor rotation rate, falling to 0.92 of the stall torque at 175 Hz. In the high-speed regime above 175 Hz, the torque falls steeply to zero at 350 Hz. The torque–speed curve is thus expressed as a piecewise linear function of the motor rotation rate, ${\mathrm{\Omega }}_m$:

$$\tau =\left\{\begin{array}{cc}{\displaystyle \left(-0.59\left(\frac{{\mathrm{\Omega }}_m}{2\pi }\right)+1300\right)}\mathrm{pN}·\mathrm{nm}\hfill & \mathrm{for}{\displaystyle 0\le \frac{{\mathrm{\Omega }}_m}{2\pi }\le 175}\mathrm{Hz}\hfill \\ {\displaystyle \left(-6.83\left(\frac{{\mathrm{\Omega }}_m}{2\pi }\right)+2392\right)}\mathrm{pN}·\mathrm{nm}\hfill & \mathrm{for}{\displaystyle 175\le \frac{{\mathrm{\Omega }}_m}{2\pi }\le 300}\mathrm{Hz}\hfill \end{array}\right.$$

(6)

Figure 3 shows the torque–speed curve as a solid black line. In each of our simulations, we ensured that the prescribed motor speed and the computed torque load formed a pair that corresponded to a point on that line.

2.3. Dynamically Similar Experiments

Experiments were performed in a 45-liter tank (0.3 m × 0.5 m × 0.5 m high) filled with incompressible silicone oil (Clearco^®) with density 970 kg/m³ and dynamic viscosity $\mu$ = $1.13\times {10}^2$ kg/(m·s) at 22 °C, about ${10}^5$ times that of water. The length and speed scales in the experiment ensured that the incompressible Stokes equations Equation (2) were valid. The viscosity of the oil drifted from the manufacturer’s stated value ($\mu$ = $1.00\times {10}^2$ kg/(m·s)) very slowly over a two-year period, so we determined the modified viscosity by measuring the torque on rotating cylinders at the center of the tank and recorded data within two months of that measurement.

The theoretical value for torque per unit length on an infinite rotating cylinder in Stokes flow is $\sigma =4\pi \mu \mathrm{\Omega }{r}^2$, where $\mu$ is the dynamic viscosity of the fluid, $\mathrm{\Omega }$ is the angular rotation rate, and r is the cylindrical radius. We measured the torque $\tau$ on a rotating cylinder with radius $r=6.35\pm 0.2$ mm and length $\ell =149\pm 1$ mm and, by assuming $\tau =\ell \sigma$, used the data to solve for the viscosity of the fluid. We also assumed that the finite size of the tank did not affect the torque value in the middle, which was more than $25r$ from the nearest boundary. Before each data collection run, we measured the temperature of the oil with a NIST-traceable calibrated thermistor (Cole-Parmer Digi-Sense-AO-37804-04 Calibrated Digital Thermometer) and adjusted the previously determined viscosity using the manufacturer’s temperature coefficient of viscosity $1.00\times {10}^{-6}$ kg/(m·s)/°C. See Section 2.3.2 for a detailed description of the torque measurements.

2.3.1. Fabricating Helices

We fabricated helices of varying wavelengths ($2.26<\lambda /R<11.88$) by wrapping straight stainless steel welding wire around cylindrical aluminum mandrels with different helical V-grooves precisely machined using a CNC lathe. The V-grooves transition to a flat face with a straight groove, to which the remaining straight section can be clamped; see Figure 4.

Mandrels were held on a lathe, and the wire was hand-spun into the V-groove. The straight sections were secured to the flat faces, which left straight stems aligned with the axes of the helices to be attached to the motor via a rigid shaft adapter. Residual tension in the wires caused the wavelengths and radii to vary after they were removed from the mandrels. The helices were forced onto a precision stainless steel rod with radius $R=6.350\pm 0.013$ mm for annealing. The helices on the rod were then placed into a tube furnace (MTI GSL-1500X) and annealed at 900 degrees Celsius for two hours, which removed most of the variation in the radii of the helices and fixed the helical wavelengths. The helix parameters used in the experiments are listed

Table 3. Wavelengths and lengths of helices.

$\mathit{\lambda }/\mathit{R}$	$\mathit{L}/\mathit{R}$
2.26 ± 0.13	22.3 ± 0.5
3.88 ± 0.01	24.3 ± 0.5
5.86 ± 0.08	30.0 ± 0.5
8.65 ± 0.01	23.3 ± 0.5
10.91 ± 0.01	24.2 ± 0.5
11.88 ± 0.01	23.1 ± 0.5

The helical wavelength $\lambda$ and axial length L are expressed in terms of the helical radius R; $R=6.35\pm 0.10$ mm in all cases. The filament radius was $a/R=0.111$ for all helices.

:

2.3.2. Axial Torque Measurements

To measure the dependence of torque on boundary distance, we secured the tank onto a horizontal stage that allowed for motion in the x-direction, as shown in Figure 5. The motion of the stage was controlled by a linear guide with a worm gear screw that advanced the stage 0.3 mm per revolution. The screw was turned using a computer-controlled NEMA 23 stepper motor with a resolution of 400 steps/rev. This gave better than 100 μm precision in controlling the boundary distance, which was necessary: the step size near the boundary was as small as 0.5 mm.

Torque measurements were made for both cylinders and helices using similar methods. The objects were held in a rigid shaft adapter and then lowered until centered in the tank using a vertical translation stage built from 80-20^® extruded aluminum.

At the beginning of each data set, we first adjusted the vertical tilt of the object until it was parallel to the boundary. Next, we manually adjusted the horizontal stage so that the cylinder or helix touched the front vertical boundary of the tank. We used total internal reflection to form an image of the object that could be used as a reference to find where the edge of the object just contacted the boundary, which occurs when the image appears to touch the object.

The torque was measured using a FUTEK TFF400, 10 in-oz, Reaction Torque Sensor. The cylinder and helices were driven by a variable speed DC motor with a magnetic encoder (Pololu 298:1 Micro Metal Gear Motor with Magnetic Encoder) and housed inside of a 3D-printed enclosure that included sleeve bearings to minimize frictional torque. The power and signal wires were fed through a 6.32 mm opening at the center of the torque sensor. The wires were then fixed to the outside structure so that they did not create a torque when measurements were taken. The encoder output was read by the counter input on a National Instruments USB6211 M series multifunction DAQ. The torque signal was amplified using an amplifier/driver (Omega DP25B-E-A 1/8 DIN Process Meter and Controller) and its output fed into the same National Instruments data acquisition board’s analog to digital input with a resolution of 250 thousand samples per second, which is much faster than any time scales in the experiment.

Data were taken with the DC motor rotating at varying speeds and with the objects located at a distance from the boundary set by the horizontal stage. The torque and motor frequency were simultaneously recorded using MATLAB to acquire and plot them. We used MATLAB and a motor controller (ARDUINO MEGA 2560 with an ADAFRUIT Motor Shield v.2) to control the motor. However, the motor rotation varied depending on the axial load, so we divided the signal from the torque sensor by the frequency data from the counter input to get the torque per unit frequency at each boundary distance, see Figure 6. A MATLAB data acquisition GUI included the temperature and distance values, ensuring that the acquisition parameters were stored with the raw data.

Data were taken for approximately 60 rotation periods for both CW and CCW rotation at each boundary location. The frequency signal occasionally showed large spikes that affected the average torque-per-frequency value because the torque signal did not show a corresponding jump. We considered this to be the result of the encoder miscounting the rotation rate or the counter input in the DAQ misreading the signal from the encoder. We used MATLAB’s outliers function to remove such frequency spikes that were more than nine median absolute deviations from the median calculated in a moving window ten data points wide and replaced them with the average of the adjoining data points. The number of outliers was less than 1% of the data points, so this frequency smoothing should not have biased the averaging significantly.

The difference between mean CW and CCW rotation values, which should have been the same, was used to establish the uncertainty in the experimental measurements. An analysis script read the geometric parameters and data files for a given set of measurements (cylinder or helix) and plotted the data versus boundary distance. We scaled the torque using a unit of $[\mu \mathrm{\Omega }{r}^2\ell ]$ (cylinder) or $\left[\mu \mathrm{\Omega }{R}^{2}L\right]$ (helix), where $\mu$ is the fluid viscosity, $\mathrm{\Omega }$ is the angular speed, r is the cylindrical radius, ℓ is the cylindrical length, R is the helical radius, and L is the helical axial length. Plots of the dimensionless torque for cylinders and helices are shown in Figure 7, Figure 8 and Figure 9. Using these units for the torque allowed for easy comparison between experiments, numerical simulations, and theory.

2.4. Summary of Algorithms and Data Analysis

Two separate sets of simulations are presented in this paper. For those with a helix model or a bacterium model, the results were averaged over 16 evenly spaced phases as described in Equation (1) of the flagellar centerline.

(i) The goal of the first set of simulations was to calibrate the MRS and MIRS methods by finding the optimal factors (${\gamma }_c$ for a cylindrical cell body and ${\gamma }_f$ for a helical flagellum) and the optimal regularization parameters (${ϵ}_c$ and ${ϵ}_f$), as reported in Table 2. Equation (3) was used to solve for the force ${\mathbf{f}}_k$ at each discretized point ${\mathbf{x}}_k$ in a free space, whereas Equation (4) was used for simulations near a plane wall. The resulting net torque of each rotating structure was then compared with the results from theory for a cylinder or from experiments for a helix, as described in Section 3.1.

(ii) The goal of the second set of simulations was to assess the motility performance of the force-free and torque-free bacterium models with boundary effects incorporated.

Step 1: Equation (5) was used with $S_{ϵ}$ (for simulations in a free space) or with $S_{ϵ}^{*}$ (for simulations with a plane wall). Different combinations of the cell body size, flagellar wavelength, and distance to the wall were simulated. We used five values for the length ℓ and five values for the radius r shown in Table 2. These values are within the range of normal E. coli [21]. We used 18 wavelengths $\lambda$ that cover a range of biological values (2.22 ± 0.2 μm) and values that are shorter and longer than the biological values (Table 2 and Figure 2). The set of geometric parameters, together with 22 distance values d measured from the flagellar axis of symmetry to the wall, resulted in 9900 simulations. From each simulation, we obtained the axial component of the translational velocity U, the magnitude of the axial-component of the hydrodynamic drag on the cell body F, and the magnitude of the axial-component of the hydrodynamics torque on the cell body $\tau$. For each body geometry (450 total), we performed a simulation in free-space to ensure the convergence of MIRS calculations to MRS calculations as the distance $d\to \infty$.

Step 2: The torque value $\tau$ was output from each simulation in Step 1 with the motor frequency set to $154\mathrm{Hz}$. That torque-frequency pair was then used to determine the load line and its intersection with the torque–speed, as discussed in Section 2.2 and shown in Figure 3. Each motor frequency ${\mathrm{\Omega }}_m/2\pi$ on the torque–speed curve was given as some multiple q of $154\mathrm{Hz}$. The simulation outputs were scaled by q, since they were all linear with motor frequency; i.e., $\left(U,F,\tau \right)\to q\left(U,F,\tau \right)$. These scaled quantities were then used to calculate the performance measures. Results are presented in Section 3.2 and Section 3.3.

3. Results

3.1. Verifying the Numerical Model and Determining the Optimal Regularization Parameters

When using MRS or MIRS, the choice of the regularization parameter for a given discretization (cylinder) or filament radius (helix) of the immersed structure has generally been made without precise connection to real-world experiments, because there are large uncertainties in biological and other small-scale measurements. We therefore used theory, as described below, and dynamically similar experiments, as described in Section 2.3, to determine the optimal regularization parameters for the two geometries used in our bacterial model: a cylinder and a helix.

3.1.1. Finding the Optimal Regularization Parameter for a Rotating Cylinder

Jeffrey and Onishi (1981) derived a theory for the torque per length on an infinite cylinder rotating near an infinite plane wall [27] that was used previously to calibrate numerical simulations of helical flagella [24]. The torque per unit length $\sigma$ on an infinite cylinder is given as:

$$\sigma =4\pi \mu \mathrm{\Omega }{r}^2\frac{d}{{(d^2-r^2)}^{1/2}}$$

(7)

where $\mu$ is the dynamic viscosity of the fluid, $\mathrm{\Omega }$ is the angular rotation speed, r is the cylindrical radius, and d is the distance from the axis of symmetry to the plane wall.

We used this theoretical value as a common reference point between the experiments and simulations to establish optimal computational parameters, but note that this theory has not been experimentally tested outside of the present work. We assumed Equation (7) is valid for our experiments and simulations, though this assumption as applied to experiments ignored the finite size of the tank. To control for end effects in the experiments, we measured the torque with only the first 3 cm inserted into the fluid and with the full cylinder inserted at the same boundary locations. We subtracted the torque found for the short section from the torque found for the full insertion of the cylinder. In simulations, we controlled for finite-length effects by measuring the torque on a middle subsection of the simulated cylinder, as discussed below.

Our experimental data are shown in Figure 7, with the torque made dimensionless using the quantity $\mu \mathrm{\Omega }{r}^2\ell$, where $\mu$ is the fluid viscosity, $\mathrm{\Omega }$ is the rotation rate, r is the cylindrical radius, and ℓ is the cylindrical length. The mean squared error (MSE) between experiments and theory is MSE $\le 6\%$ when calculated for the boundary distances where $d/r>1.1$ (i.e., the distance from the boundary to the edge of the flagellum is ≥1 mm). The theory asymptotically approaches infinity as the boundary distance approaches $d/r=1$, which skewed the MSE unrealistically. For the data where $d/r\ge 2$, the mean squared error is less than 1%.

In numerical simulations of the cylinder, the computed torque value depended on both the discretization and regularization parameter. Having found good correspondence with the experiments, we used Equation (7) to find an optimal regularization parameter for a given discretization of the cylinder (see Table 2: cylinder part). The discretization size of the cylindrical model $d{s}_c$ was varied among 0.192 μm, 0.144 μm, and 0.096 μm. For each ds_c, an optimal discretization factor γ_c was found by minimizing the MSE between the numerical simulations and the theoretical values using the computed torque in the middle two-thirds of the cylinder to avoid end effects. The optimal factor was found to be ${\gamma }_c=$ 6.4 for all the discretization sizes. We used the finest discretization size for our model bacterium as reported in Table 2 since it returned the smallest MSE value of 0.36%.

3.1.2. Finding the Optimal Regularization Parameter for a Rotating Helix Far from a Boundary

Simulated helical torque values also depend on the discretization and regularization parameter, but there is no theory for a helix to provide a reference. Other researchers have determined the regularization parameter using complementary numerical simulations, but the reference simulations also have free parameters that may have affected their results [25].

Thus, we used dynamically similar experiments, as described in Section 2.3, to determine the optimal filament factor, ${\gamma }_f=2.139$, for a helix filament radius $a/R=0.111$. Torque was measured for the six helical wavelengths given in Table 3 when the helix was far from the boundary. The optimal filament factor ${\gamma }_f=2.139$ was found by the following steps: (i) varying ${ϵ}_f$ for each helix until the percent difference between the experiment and simulation was under 5%; and (ii) averaging the ${ϵ}_f$ values found in Step (i). In these simulations, the regularization parameter and discretization size are both equal to ${\gamma }_{f}a$. The results are shown in Figure 8, with the torque values non-dimensionalized by the value $\mu \mathrm{\Omega }{R}^{2}L,$ where $\mu$ is the fluid viscosity, $\mathrm{\Omega }$ is the rotation rate, R is the helical radius, and L is the axial length. The optimized simulations returned an average percent difference of $2.4\pm 1.7$% compared to the experimental values.

We checked whether helices with different filament radii could be accurately simulated using our optimized ${\gamma }_f$ to scale the regularization parameter, i.e., (${ϵ}_f={\gamma }_{f}a$), to account for relative size of the filament, as is commonly done [29,30,31,32,33]. We computed torque values that matched the experimental values given in Rodenborn et al. (2013) [26], which used a filament radius $a/R=0.063$. The results are also presented in Figure 8. The percent difference between our MRS simulations and their data is 2.5 ± 1.3%.

Martindale et al. (2016) [25] used an MRS with a surface discretization of the flagellum to calibrate their simulation parameters, whereas our MRS used a string of regularized Stokeslets along the helical centerline to reduce the computational cost in the MIRS calculations. As a final test, we used the freely available and calibrated code for the MRS with surface discretization from Rodenborn et al. (2013) [26] to compute torque values for our $a/R=0.111$ data and for their $a/R=0.063$ data. Figure 8 shows the torque comparison of the Rodenborn et al. (2013) [26] surface discretized MRS (solid curves), our centerline distribution MRS (triangles), and experiments (circles). The percent difference between their MRS and the experiments is $3.6\pm 3.4$%. The percent difference between our MRS and simulations in Rodenborn et al. (2013) is 1.8 ± 3.7%. Thus, our MRS with a centerline distribution using the optimal filament factor ${\gamma }_f$ worked very well for another filament radius and other helical wavelengths.

3.1.3. Torque on Rotating Helices Near a Boundary

To determine how boundaries affect bacterial motility, we used our optimized value for ${\gamma }_f$ in our MIRS simulations to compute the torque as a function of boundary distance, as shown in Figure 9. The computed torque values and measured torque values also show excellent agreement at most boundary distances, except for the shortest wavelength $\lambda /R=2.26$. We note that this helix had the largest variation in wavelength, as reported in Table 3. Furthermore, the torque for short wavelengths is more sensitive to variation in wavelength as compared to variation at longer wavelengths, which likely explains the difference between simulation and experiment for this geometry, whereas for the other wavelengths, the simulated values are generally within the uncertainty in the experiments for all boundary distances.

3.2. Speed Measurements to Assess Performance

The motion of bacteria through their environment enables them to find nutrients. Indeed, it has been suggested that the purpose of bacterial motility is primarily to perform chemotaxis [4]. Living in a microscopic environment where thermal effects are significant, bacteria must be able to sample chemical concentrations faster than diffusion causes those concentrations to change [4,12], so moving faster may confer a survival advantage.

The low-speed operating regime of the bacterial motor (below 175 Hz) is thermodynamically more efficient than the high-speed regime. A simple model gives the fraction of energy lost to friction in the motor as $({\tau }_0-\tau )/{\tau }_0$, where ${\tau }_0$ is the stall torque and $\tau$ is the operating torque at a given frequency [14]. In the low-speed regime, $\tau \ge 0.92{\tau }_0$, so that the power output of the motor is greater than 92% of the power input. However, the low-speed regime may be less operationally reliable for motility; the flatness of the torque–speed curve implies that small increases in load correspond to large decreases in motor rotation rate, so the bacterium risks stalling and may be unable to restart its motor. Using our simulations, we determined the swimming speed and motor rotation rate for different bacterial geometries at different distances to a solid boundary and assessed the performance of bacterial geometries typically associated with swimming.

3.2.1. Optimal Flagellar Wavelength

We first consider the effect of different flagellar wavelengths on swimming speed and motor rotation rate, as shown in Figure 10a,b. Swimming speed and motor rotation rate are shown as heat maps for different flagellar wavelengths at different distances to the boundary. The heat maps show the median values computed among all 25 bacterial body geometries we investigated (Table 2).

The maximum of all the median swimming speed values is about $26\mathsf{\mu }{\mathrm{ms}}^{-1}$, and it occurs far from the boundary for a wavelength near $8R$. For long and short flagellar wavelengths, the swimming speed at all distances is much lower than the maximum. Long wavelengths yield about $10\mathsf{\mu }{\mathrm{ms}}^{-1}$, whereas very short wavelengths give values closer to $1\mathsf{\mu }{\mathrm{ms}}^{-1}$. For the flagellar wavelength of $\lambda /R=11.1$ that is typical for E. coli, the swimming speed is about $25\mathsf{\mu }{\mathrm{ms}}^{-1}$ far from the boundary, whereas it drops to about $20\mathsf{\mu }{\mathrm{ms}}^{-1}$ very near the boundary.

Interestingly, the flagellar wavelengths that correspond to swimming speeds near the maximum in Figure 10a also correspond to motor rotation rates in the low end of the high-speed regime in the torque–speed curve, so that the motion is both thermodynamically efficient and operationally reliable. A wavelength of $\lambda /R=8$ gives 190 Hz and 183 Hz far from and near the wall, respectively, which correspond to mechanical energy outputs of about 84% and 88%. Short and long wavelengths result in a weaker performance, but for different reasons: short wavelengths operate in the low-speed regime and thus are efficient but unreliable, whereas longer wavelengths operate farther into the high-speed region and thus are reliable but inefficient.

3.2.2. Boundary Effects

To illustrate how proximity to the boundary affects swimming speed and motor rotation rate, we show line plots in Figure 10c–f of the speed and rotation rate as functions of the flagellar wavelength both far from and near the boundary. The maximum, median, and minimum values among all bacterial body geometries are shown for each boundary distance. Comparing Figure 10c,e shows that proximity to the boundary does not appreciably alter the optimal wavelength: it remains near $8R$ for all body geometries both near and far from the boundary. However, proximity to the boundary does increase the difference in the swimming speed among different bodies at a given wavelength.

Far from the boundary, the difference between the maximum and minimum swimming speeds for the optimal flagellar wavelength is 14% of the maximum value of $28\mathsf{\mu }{\mathrm{ms}}^{-1}$; near the boundary, the difference is 34% of the maximum value of $26\mathsf{\mu }{\mathrm{ms}}^{-1}$. Figure 10d,f shows the motor rotation rate is less sensitive to the body geometry and proximity to the surface than the swimming speed. Far from the boundary, the difference between the maximum and minimum rotation rates for the optimal flagellar wavelength $8R$ is 6% of the maximum value $198\mathrm{Hz}$; near the boundary, the difference is 6% of the maximum value of $190\mathrm{Hz}$.

To further probe the effect of the cell body geometry on swimming speed and motor rotation rate, we show heat maps of the speed and rotation rate fixed at the typical E. coli wavelength $\lambda /R=11.1$ as functions of the length and radius of the cylindrical cell body. Figure 11a–d shows the results. The translational speed is optimized for short, thin cell bodies (lower left-hand corner of Figure 11a,c) both near and far from the surface. Conversely, the slowest motor rotation rates (though all higher than 175 Hz), and therefore the most thermodynamically efficient, occur for long, thick cell bodies (upper right-hand corners of Figure 11b,d). Taken together, these two results suggest that balancing the need of a bacterium to move quickly with its need to be thermodynamically efficient would yield a cell body geometry somewhere between long, thick cell bodies and short, thin cell bodies. Interestingly, the center point of the heat maps shown in Figure 11 corresponds to the mean size of the E. coli cell body.

3.3. Energy Cost Measures to Assess Performance

The energy cost required to move is another way to assess the performance of the bacterial motility system. Here, we present simulation results of three different energy cost measures. The first measure we consider is what we term the Purcell inefficiency ${\mathcal{E}}_{Purcell}^{-1}$ given by,

$${\mathcal{E}}_{Purcell}^{-1}=\frac{\tau {\mathrm{\Omega }}_m}{FU},$$

(8)

where $\tau$ is the motor torque (or the torque on the cell body or the flagellum), ${\mathrm{\Omega }}_m$ is the motor rotation rate, F is the drag force on the cell body (or on the flagellum), and U is the swimming speed of the bacterium.

Thus, the Purcell inefficiency measures the mechanical power ($T{\mathrm{\Omega }}_m$ ) required to swim at speed U relative to the least power ($FU$) needed to translate the cell body at speed U. The Purcell inefficiency is useful because, under certain simplifying assumptions [34], it can be expressed as a function of the geometry of the cell body and the flagellum alone. The difficulty with this measure is that it does not depend on the rotation rate of the motor because all four quantities appearing in Equation (8) scale with the motor frequency (see Equation (5)). Therefore, the Purcell inefficiency cannot assess how swimming performance depends on the torque–speed characteristics of the motor and thus omits an important element of the bacterial motility system that is subject to selective forces.

The second measure is the energy cost to travel a unit distance given by

$$E^{*}=\frac{\tau {\mathrm{\Omega }}_m}{U}.$$

(9)

Several authors [14,17] have considered the distance traveled per energy output by the motor, which is the inverse of the measure we consider here. The merit of the energy cost per distance measure is that it expresses the amount of energy used by the bacterium to perform a biologically relevant task; namely, to swim one unit distance. Another advantage is that it depends on the motor rotation rate and thus can probe the effect of the torque–speed characteristics of the motor. However, it does not account for the size of the bacterium, and thus does not measure the energy cost relative to the overall metabolic budget of the organism.

To account for the metabolic energy cost required to swim a unit distance, we introduce a third measure,

$$\mathcal{E}=\frac{\tau {\mathrm{\Omega }}_m}{mU}.$$

(10)

The mass m associated with each bacterial model is $m=1.1\times {10}^{-15}\left(\pi r^{2}l\right)\mathrm{kg}$, where r is the body radius and ℓ is the body length, both measured in $\mathsf{\mu }\mathrm{m}$. Though this energy cost measure has not been considered in the literature, it was suggested earlier by Purcell [4].

3.3.1. Optimal Wavelength

We first consider the optimal flagellar wavelength predicted by the three energy cost measures, as shown in Figure 12. The top row a-c shows heat maps of the three energy cost measures as functions of flagellar wavelength and boundary distance, which correspond to the median values computed for all body geometries listed in Table 2. All three measures give an optimal wavelength near $\lambda /R=8$ (where each energy cost measure is minimal). However, the three measures differ in other ways. The Purcell inefficiency predicts that swimming near the boundary is less inefficient than swimming far from the boundary, whereas the opposite is true for the energy per distance and metabolic cost measures. At a wavelength of $8R$, the minimum Purcell inefficiency value is about 84 (or $1/84=1.2\%$ if calculated as Purcell efficiency), the minimum energy per distance measure is $5.0\times {10}^{-11}{\mathrm{Jm}}^{-1}$, and the minimum metabolic energy cost is $3.1\times {10}^4\mathrm{J}{\mathrm{m}}^{-1}{\mathrm{kg}}^{-1}$.

3.3.2. Boundary Effects

To evaluate how proximity to the surface affects the predictions of the energy cost measures, we show line plots in Figure 12 of the measures as functions of flagellar wavelength far from ($d/R=51$) in d-f and near the boundary ($d/R=2.74$) in g–i. The maximum, median, and minimum values among all body geometries are shown for each wavelength. The Purcell inefficiency is the least sensitive of the measures to changes in the body size. For a wavelength of $\lambda /R=8$, the difference between the maximum and the minimum is 8% of the maximum (110 vs. 101). Near the boundary, the difference increases to 13% of the maximum value (94 vs. 82).

The energy per distance measure is more sensitive to the body size, and the sensitivity increases near the boundary. For a wavelength of $\lambda /R=8$, the difference between the maximum and minimum values is 16% of the maximum value of $5.5\times {10}^{-11}\mathrm{J}{\mathrm{m}}^{-1}$ far from the boundary. Near the boundary, the difference increases to 35% of the maximum value of $7.5\times {10}^{-11}\mathrm{J}{\mathrm{m}}^{-1}$. The metabolic energy cost is the measure most sensitive to the body size, though interestingly, the sensitivity decreases with proximity to the boundary. At a wavelength of $\lambda /R=8$, the difference between the maximum and minimum value far from the boundary is 51% of the maximum value of $4.5\times {10}^4\mathrm{J}{\mathrm{m}}^{-1}{\mathrm{kg}}^{-1}$. Near the boundary, the difference decreases to 38% of the maximum value of $5.0\times {10}^4\mathrm{J}{\mathrm{m}}^{-1}{\mathrm{kg}}^{-1}$.

Finally, we consider how the energy cost measures depend on body radius and body length at different distances to the boundary. In Figure 13, we show heat maps of the three energy cost measures fixed at the typical E. coli wavelength $\lambda /R=11.1$, as functions of the radius and length. The Purcell inefficiency shown in Figure 13 gives different optimal body geometries near and far from the boundary: far from the boundary, the short, thick cylinders (top left corner of Figure 13a) are the least inefficient; near the boundary, the short, thin cylinders (bottom left corner of Figure 13d) are the least inefficient. The energy per distance measure gives the same optimal body far from and near the boundary: the lowest energy per distance cost measure is given by short, thin cylinders (bottom left corners of Figure 13b,e). The metabolic cost measure gives the same optimal body near and far from the surface, though it is opposite of the optimal body predicted by the energy per distance measure: the lowest metabolic cost measure occurs for cylinders that are long and thick (top right corners of Figure 13c,f).

4. Discussion

In this work, we used the method of images for regularized Stokeslets (MIRS) to simulate a motile flagellate bacterium moving near a solid boundary. We determined the regularization parameter in the method by conducting dynamically similar macroscopic experiments with rotating cylinders and rotating helices near a solid boundary and comparing the results to equivalent simulations. By varying the regularization parameters, we were able to find optimal values that matched the experimental results within 5%. Having calibrated MIRS, we simulated various bacterial morphologies to assess their swimming performance. We assessed swimming performance using multiple measures: swimming speed, motor rotation rate, the Purcell inefficiency, energy cost per distance, and metabolic energy cost. An important and novel addition to our simulations was to ensure that the torque and rotation rate match a point on the experimentally measured torque–speed curve [18] in all of our calculated measures.

Using our MIRS calibration method, we found that the optimal discretization factor for a cylinder is ${\gamma }_c=6.4$ for the surface discretizations we used, which may be used as a reference value for other researchers who simulate rotating cylinders using MRS or MIRS. We also found an optimal filament factor ${\gamma }_f=2.139$ when using MRS and MIRS with each helix modeled as a string of regularized Stokeslets along the helix centerline. Our computations showed excellent agreement with the experimental helix torque values at most boundary distances (Figure 8 and Figure 9). Selecting an appropriate regularization parameter for a center-line discretization of helices with MIRS has been considered by other researchers. Martindale et al. (2016) [25] benchmarked their center-line discretization of a helix with a surface discretization model. They reported that the optimal filament factor should be in the range $1\le {\gamma }_f\le 3$ to keep the percent difference less than about 10% in their simulations, which is consistent with our results.

In MRS/MIRS, using a centerline distribution for a model helix (or flagellum) with a calibrated regularization parameter is more useful than a surface discretization for several reasons: (i) the computational cost is significantly reduced because the matrix system for the centerline distribution is much smaller than for a surface discretization; (ii) simulations of very short helical wavelengths using a centerline distribution do not encounter discretization issues such as overlapping cross-sections; (iii) in a centerline distribution, the point connecting the cylindrical cell body and the tapered helical flagellum can be considered as the motor location, whereas the motor location in a surface discretization is hard to define because of the small gap between the cell body and the flagellum needed to allow counter-rotation between the cell body and the flagellum.

Interestingly, all five performance measures we computed with our calibrated model—swimming speed, motor speed, Purcell inefficiency, energy per distance, and metabolic energy cost—predict an optimal flagellar wavelength of $\lambda /R\approx 8$, where R is the helical radius of the flagellum. This result agrees with the work of Zhang et al. (2014) [10] who studied the Purcell efficiency of a rotating helix, whereas our model includes a cell body with rotation and translation. Furthermore, this prediction occurs both near and far from the surface and for all body geometries, which suggests that the bacterial wavelengths may be selected independently of body shape or surface proximity. Although the five measures agreed on the optimal flagellar wavelength, they disagreed on the optimal body size.

We can propose an experiment to determine which performance measures best predict which bacterial attributes are most efficient. From an evolutionary standpoint, the most commonly encountered traits in a population should correspond to those traits that confer optimal energetic efficiency in bacterial motility. A comparison of the distribution of traits in a representative population to the traits predicted by the energy cost measures to be most efficient would reveal which of the measures has the greatest predictive power.

Among the efficiency measures we studied, the Purcell efficiency has been the most popular quantity of analysis, but we believe it has several important shortcomings that warrant discussion. First, the Purcell efficiency is dependent only on the geometry of the body and flagellum and not on the motor’s torque–speed response characteristics. From a physical standpoint, it is interesting to find such an invariant quantity, but from a biological standpoint, it does not assess the bacterial motility system’s thermodynamic efficiency because it ignores motor mechanics.

Second, the Purcell efficiency is defined to be the ratio between the minimum power required to translate the cell body and the power actually dissipated during the bacterial motion. In our simulations, we find the maximal efficiency is in the range of 1–2%, similar to what others have found [4,10,13]. These two quantities (the minimum power vs. the actual power) are clearly of very different orders, which implies that the least power needed may not be an appropriate reference quantity. To provide a biophysical interpretation to the least power needed to translate the cell body, some authors have suggested that it represents the “useful” portion of the power dissipated during motion [9,13], but we believe this is a misconception. The bacterium is non-inertial; therefore, the force acting on the cell body by the fluid is exactly balanced by the force acting on the flagellum by the fluid (assuming no net body forces). Both the bacterial body and the flagellum have the same axial velocity (in a rigid model); therefore, the power dissipated due to the axial fluid drag on the body is exactly compensated by the power input by the axial fluid force exerted on the flagellum.

Finally, as Purcell noted in 1977, the efficiency of the bacterial motility system is probably best characterized by the energy consumption relative to the overall metabolic budget of the organism [4]. This suggestion led us to consider the metabolic energy cost $\mathcal{E}$ (Equation (10)) introduced in this paper. The actual amount of that metabolic budget used for motility is a small fraction, so Purcell [4] also suggested that bacterial motility is not really subject to strong selective forces toward optimal efficiency. Our data do not say whether evolutionary processes tend to minimize the energy cost of bacterial motility, but a plausible counterargument is that the bacterium needs to consume most of its energy for other biological functions and has only a small fraction available for motility. Thus, small absolute changes in energy consumption correspond to large relative changes in the energy available for motility, resulting in a significant selective pressure to make the motility system as efficient as possible.

Many research questions about how physical interactions between bacteria and their environment result in selective pressures in evolutionary processes remain open, despite significant progress in the field. Modern computational simulations and methods such as MRS and MIRS will remain important tools for quantifying microscopic bacterial motion with precision. Calibrating these methods with experimental results in cases where exact theory is unavailable can help to ensure simulations give accurate quantitative results. In future work, we will extend the macroscopic experimental system to consider a wider variety of possible geometries relevant to bacterial motility and make comparisons with biological measurements.