The N-flagella problem: Elastohydrodynamic motility transition of multi-flagellated bacteria

TL;DR

This paper presents a theoretical analysis of how multi-flagellated bacteria transition from non-motile to swimming states, revealing that the most unstable mode leads to linear swimming combined with rotation, consistent with experimental observations.

Contribution

The study analytically solves the elastohydrodynamic stability problem for multi-flagellated bacteria, identifying the dominant motility mode and its relation to bacterial morphology.

Findings

Most unstable mode causes linear swimming with rotation

Up to 6 unstable modes related to rigid-body motion

Results align with experimental bacterial swimming behavior

Abstract

Peritrichous bacteria such as Escherichia coli swim in viscous fluids by forming a helical bundle of flagellar filaments. The filaments are spatially distributed around the cell body to which they are connected via a flexible hook. To understand how the swimming direction of the cell is determined, we theoretically investigate the elastohydrodynamic motility problem of a multi-flagellated bacterium. Specifically, we consider a spherical cell body with a number N of flagella which are initially symmetrically arranged in a plane in order to provide an equilibrium state. We analytically solve the linear stability problem and find that at most 6 modes can be unstable and that these correspond to the degrees of freedom for the rigid-body motion of the cell body. Although there exists a rotation-dominated mode that generates negligible locomotion, we show that for the typical morphological parameters of bacteria the most unstable mode results in linear swimming in one direction accompanied by rotation around the same axis, as observed experimentally.

Equations

$$d{\mathbf{F}}^{(i)}(s)=\left[C_t{\mathbf{t}}^{(i)}{\mathbf{t}}^{(i)}+C_n(1-{\mathbf{t}}^{(i)}{\mathbf{t}}^{(i)})\right]\cdot {\mathbf{u}}^{(i)}(s)ds,$$

$${\mathbf{u}}^{(i)}=\mathbf{U}+\mathbf{\Omega }\times ({\mathbf{x}}^{(i)}+{\mathbf{\xi }}^{(i)})+{\mathbf{\omega }}^{(i)}\times {\mathbf{\xi }}^{(i)},$$

$${\mathbf{\xi }}^{(i)}(s)={\int }_0^s{\mathbf{t}}^{(i)}(s^{\prime })d{s}^{\prime }$$

$${\mathbf{F}}^{(i)}=\left[{\int }_0^{L^{(i)}}{\mathsf{C}}_{i}ds\right]\cdot \mathbf{U}+\left[{\int }_0^{L^{(i)}}{\mathsf{C}}_i\cdot {\mathsf{A}}_{i}ds\right]\cdot \mathbf{\Omega }+\left[{\int }_0^{L^{(i)}}{\mathsf{C}}_i\cdot {\tilde{\mathsf{A}}}_{i}ds\right]\cdot {\mathbf{\omega }}^{(i)}.$$

$${\left[{\mathsf{A}}_i\right]}_{pq}={\epsilon }_{pqr}(x_r^{(i)}+{\xi }_r^{(i)})\text{and}{\left[{\tilde{\mathsf{A}}}_i\right]}_{pq}={\epsilon }_{pqr}{\xi }_r^{(i)},$$

$${\mathbf{F}}_{body}+\sum _{i=1}^N{\mathbf{F}}^{(i)}=0.$$

$${\mathbf{M}}^{(i)}=\left[{\int }_0^{L^{(i)}}{\mathsf{A}}_i^T\cdot {\mathsf{C}}_{i}ds\right]\cdot \mathbf{U}+\left[{\int }_0^{L^{(i)}}{\mathsf{A}}_i^T\cdot {\mathsf{C}}_i\cdot {\mathsf{A}}_{i}ds\right]\cdot \mathbf{\Omega }+\left[{\int }_0^{L^{(i)}}{\mathsf{A}}_i^T\cdot {\mathsf{C}}_i\cdot {\tilde{\mathsf{A}}}_{i}ds\right]\cdot {\mathbf{\omega }}^{(i)},$$

$${\mathbf{M}}_{body}+\sum _{i=1}^N{\mathbf{M}}^{(i)}=0.$$

$${\tilde{\mathbf{M}}}^{(i)}=\left[{\int }_0^{L^{(i)}}{\tilde{\mathsf{A}}}_i^T\cdot {\mathsf{C}}_{i}ds\right]\cdot \mathbf{U}+\left[{\int }_0^{L^{(i)}}{\tilde{\mathsf{A}}}_i^T\cdot {\mathsf{C}}_i\cdot {\mathsf{A}}_{i}ds\right]\cdot \mathbf{\Omega }+\left[{\int }_0^{L^{(i)}}{\tilde{\mathsf{A}}}_i^T\cdot {\mathsf{C}}_i\cdot {\tilde{\mathsf{A}}}_{i}ds\right]\cdot {\mathbf{\omega }}^{(i)}.$$

$${\mathbf{M}}_{elast}^{(i)}=-{\kappa }^{(i)}{\theta }^{(i)}{\mathbf{e}}_{\perp }^{(i)},$$

$${\mathsf{P}}_i\cdot {\tilde{\mathbf{M}}}^{(i)}+{\mathbf{M}}_{elast}^{(i)}=0,$$

$$\mathcal{A}\left(\begin{array}{c}\mathbf{U}\\ \mathbf{\Omega }\\ {\mathbf{\omega }}^{(i)}\end{array}\right)=\left(\begin{array}{c}0\\ 0\\ -{\mathbf{M}}_{elast}^{(i)}\end{array}\right),$$

$$\mathcal{A}=\left(\begin{array}{ccc}{C}_{D}1+\sum _{i=1}^N{\int }_0^{L^{(i)}}{\mathsf{C}}_{i}ds& \sum _{i=1}^N{\int }_0^{L^{(i)}}{\mathsf{C}}_i\cdot {\mathsf{A}}_{i}ds& {\int }_0^{L^{(i)}}{\mathsf{C}}_i\cdot {\tilde{\mathsf{A}}}_{i}ds\\ \sum _{i=1}^N{\int }_0^{L^{(i)}}{\mathsf{A}}_i^T\cdot {\mathsf{C}}_{i}ds& C_{R}1+\sum _{i=1}^N{\int }_0^{L^{(i)}}{\mathsf{A}}_i^T\cdot {\mathsf{C}}_i\cdot {\mathsf{A}}_{i}ds& {\int }_0^{L^{(i)}}{\mathsf{A}}_i^T\cdot {\mathsf{C}}_i\cdot {\tilde{\mathsf{A}}}_{i}ds\\ {\mathsf{P}}_i\cdot {\int }_0^{L^{(i)}}{\tilde{\mathsf{A}}}_i^T\cdot {\mathsf{C}}_{i}ds& {\mathsf{P}}_i\cdot {\int }_0^{L^{(i)}}{\tilde{\mathsf{A}}}_i^T\cdot {\mathsf{C}}_i\cdot {\mathsf{A}}_{i}ds& {\mathsf{P}}_i\cdot {\int }_0^{L^{(i)}}{\tilde{\mathsf{A}}}_i^T\cdot {\mathsf{C}}_i\cdot {\tilde{\mathsf{A}}}_{i}ds\end{array}\right).$$

$$\left(\begin{array}{cccc}{\mathsf{K}}_{TT}& {\mathsf{K}}_{TR}& {\mathsf{K}}_{TF}^{(i)}& \\ {\mathsf{K}}_{RT}& {\mathsf{K}}_{RR}& {\mathsf{K}}_{RF}^{(i)}& \\ {\mathsf{P}}_i\cdot {\mathsf{K}}_{FT}^{(i)}& {\mathsf{P}}_i\cdot {\mathsf{K}}_{FR}^{(i)}& {\mathsf{P}}_i\cdot {\mathsf{K}}_{FF}^{(i)}& \end{array}\right)\left(\begin{array}{c}\mathbf{U}\\ \mathbf{\Omega }\\ {\mathbf{\omega }}_n^{(i)}\end{array}\right)=\left(\begin{array}{c}-{\sum }_{i=1}^N{\mathsf{K}}_{TF}^{(i)}\cdot {\mathbf{\omega }}_t^{(i)}\\ -{\sum }_{i=1}^N{\mathsf{K}}_{RF}^{(i)}\cdot {\mathbf{\omega }}_t^{(i)}\\ -{\mathbf{M}}_{elast}^{(i)}-{\mathsf{P}}_i\cdot {\mathsf{K}}_{FF}^{(i)}\cdot {\mathbf{\omega }}_t^{(i)}\end{array}\right).$$

$$\mathsf{K}={\mathsf{R}}_i\cdot {\mathsf{K}}_i\cdot {\mathsf{R}}_i^{-1}.$$

$$\mathsf{K}=\mathsf{R}\left(-\alpha ;{\mathbf{e}}^{(i)}\right)\cdot \mathsf{K}\cdot \mathsf{R}\left(\alpha ;{\mathbf{e}}^{(i)}\right),$$

$${\mathsf{K}}_i=\left(\begin{array}{ccc}{K}_{11}& K_{12}& 0\\ K_{21}& K_{22}& 0\\ 0& 0& K_{33}\end{array}\right),$$

$$\left(\begin{array}{cc}{K}_{11}& K_{12}\\ K_{21}& K_{22}\end{array}\right)=\left(\begin{array}{cc}{K}_s& 0\\ 0& K_s\end{array}\right)+\left(\begin{array}{cc}0& K_{ss}\\ -K_{ss}& 0\end{array}\right),$$

$${\mathsf{P}}_i\cdot \mathsf{K}=\mathsf{K}\cdot {\mathsf{P}}_i.$$

$$\left(\begin{array}{ccc}{\mathsf{K}}_{TT}& {\mathsf{K}}_{TR}& {\mathsf{K}}_{TF}^{(i)}\\ {\mathsf{K}}_{RT}& {\mathsf{K}}_{RR}& {\mathsf{K}}_{RF}^{(i)}\\ {\mathsf{P}}_i\cdot {\mathsf{K}}_{FT}^{(i)}& {\mathsf{P}}_i\cdot {\mathsf{K}}_{FR}^{(i)}& {\mathsf{K}}_{FF}^{(i)}\end{array}\right)\left(\begin{array}{c}\mathbf{U}\\ \mathbf{\Omega }\\ {\mathbf{\omega }}_n^{(i)}\end{array}\right)=\left(\begin{array}{c}-{\sum }_{i=1}^N{\mathsf{K}}_{TF}^{(i)}\cdot {\mathbf{\omega }}_t^{(i)}\\ -{\sum }_{i=1}^N{\mathsf{K}}_{RF}^{(i)}\cdot {\mathbf{\omega }}_t^{(i)}\\ -{\mathbf{M}}_{elast}^{(i)}\end{array}\right).$$

$$\left(\begin{array}{ccc}{\mathsf{K}}_{TT}& {\mathsf{K}}_{TR}& {\mathsf{R}}_i\cdot {\mathsf{K}}_{TF}^{(0)}\\ {\mathsf{K}}_{RT}& {\mathsf{K}}_{RR}& {\mathsf{K}}_{RF}^{(i)}\cdot {\mathsf{R}}_i\\ {\mathsf{R}}_i\cdot {\mathsf{P}}_0\cdot {\mathsf{K}}_{FT}^{(0)}\cdot {\mathsf{R}}_i^{-1}& {\mathsf{P}}_i\cdot {\mathsf{K}}_{FR}^{(i)}& {\mathsf{R}}_i\cdot {\mathsf{K}}_{FF}^{(0)}\end{array}\right)\left(\begin{array}{c}\mathbf{U}\\ \mathbf{\Omega }\\ {\tilde{\mathbf{\omega }}}_n^{(i)}\end{array}\right)=\left(\begin{array}{c}-{\sum }_{i=1}^N{\mathsf{R}}_i\cdot {\mathsf{K}}_{TF}^{(0)}\cdot {\tilde{\mathbf{\omega }}}_t^{(i)}\\ -{\sum }_{i=1}^N\left({\mathsf{K}}_{RF}^{{}^{\prime }(i)}\cdot {\mathsf{R}}_i+{\mathsf{R}}_i\cdot {\mathsf{K}}_{FF}^{(0)}\right)\cdot {\tilde{\mathbf{\omega }}}_t^{(i)}\\ -{\mathsf{R}}_i\cdot {\tilde{\mathbf{M}}}_{elast}^{(i)}\end{array}\right),$$

$$\left(\begin{array}{ccc}{\mathsf{K}}_{TT}& {\mathsf{K}}_{TR}& {\mathsf{R}}_i\cdot {\mathsf{K}}_{TF}^{(0)}\\ {\mathsf{K}}_{RT}& {\mathsf{K}}_{RR}& {\mathsf{K}}_{RF}^{(i)}\cdot {\mathsf{R}}_i\\ {\mathsf{P}}_0\cdot {\mathsf{K}}_{FT}^{(0)}\cdot {\mathsf{R}}_i^{-1}& {\mathsf{P}}_0\cdot {\mathsf{R}}_i^{-1}\cdot {\mathsf{K}}_{FR}^{(i)}& {\mathsf{K}}_{FF}^{(0)}\end{array}\right)\left(\begin{array}{c}\mathbf{U}\\ \mathbf{\Omega }\\ {\tilde{\mathbf{\omega }}}_n^{(i)}\end{array}\right)=\left(\begin{array}{c}-{\mathbf{F}}_{prop}\\ -{\mathbf{M}}_{prop}\\ -{\tilde{\mathbf{M}}}_{elast}^{(i)}\end{array}\right).$$

$${\mathbf{F}}_{prop}$$

$${\mathbf{M}}_{prop}$$

$${\mathsf{K}}_C^{(0)}=\left(\begin{array}{ccc}{k}_C& 0& 0\\ 0& k_C& 0\\ 0& 0& K_C\end{array}\right),{\mathsf{K}}_{TF}^{(0)}=\left(\begin{array}{ccc}{k}_D& k_T& 0\\ -k_T& k_D& 0\\ 0& 0& K_T\end{array}\right),{\mathsf{K}}_{FF}^{(0)}=\left(\begin{array}{ccc}{k}_F& 0& 0\\ 0& k_F& 0\\ 0& 0& K_F\end{array}\right),$$

$${\mathsf{R}}_i$$

$${\mathbf{F}}_{prop}=\left(\sum _{i=1}^N{\omega }_0{\mathsf{R}}_i\right)\cdot {\mathsf{K}}_{TF}^{(0)}\cdot {\mathbf{e}}_z={\omega }_0{K}_T\left(-\sum _{i=1}^N{\theta }_i\sin {\Theta }_i{\mathbf{e}}_x+\sum _{i=1}^N{\theta }_i\cos {\Theta }_i{\mathbf{e}}_y+\sum _{i=1}^N{\varphi }_i{\mathbf{e}}_z\right),$$

$${\mathbf{M}}_{prop}=$$

$$+{\omega }_0{K}_F\left(-\sum _{i=1}^N{\theta }_i\sin {\Theta }_i{\mathbf{e}}_x+\sum _{i=1}^N{\theta }_i\cos {\Theta }_i{\mathbf{e}}_y+\sum _{i=1}^N{\varphi }_i{\mathbf{e}}_z\right).$$

$$\mathcal{A}\mathbf{\Phi }=\mathbf{b},$$

$$\mathbf{b}=(-{\mathbf{F}}_{prop}^T,-{\mathbf{M}}_{prop}^T,\kappa {\theta }_1,\cdots ,\kappa {\theta }_N,\kappa {\varphi }_1,\cdots ,\kappa {\varphi }_N{)}^T.$$

Full text

††preprint: APS/123-QED

I Introduction

Bacteria, which constitute the largest domain of prokaryotic microorganisms, have survived for billions of years due to their sophisticated structures lighthill1976 ; bray2000 ; phillips2012 . Many bacteria are motile and various forms of bacterial motility have been discovered, including gliding, twitching and swarming youderian1998 ; desiel2001 ; kearns2010 ; bergweb . Above all, the most common form is swimming and both its hydrodynamics basis and the chemotactic behaviour of cells have been investigated for many decades purcell1977 ; berg2003 ; berg2004 ; tindall2012 ; lauga2016 .

Bacterial swimming is achieved using propulsion from helical appendages, called flagellar filaments, attached to the cell body (typical diameter $\approx 1-2~{}$µm). A flagellar filament is a slender polymer, made of a single protein called flagellin, which takes the form of a rigid helix (typical length $\approx 10~{}$µm and diameter $\approx 40$ nm), driven in rotation by a bacterial rotary motor located at the base end of the filament. The numbers and positions of the flagellar filaments can vary greatly from cell to cell, but so-called peritrichous bacteria possess multiple flagella effectively randomly distributed on the cell surface. This group of bacteria includes well-studied organisms such as Bacillus subtilis, Salmonella enterica and the most popular model bacterium, Escherchia coli berg2004 ; guttenplan2013 ; lauga2016 .

The behaviour of a rotary motor is regulated by inter-cellular signalling proteins. When a peritrichous bacterium is swimming (so-called “run”), the distributed flagellar filaments rotate in the same direction and gather together in a helical bundle, generating essentially linear propulsion. When at least one of the rotary motors counter-rotates, the bundle of flagellar filaments comes apart and the cell changes its orientation macnab1977 (“tumble”). The resulting well-studied run-and-tumble mechanism allows peritrichous bacteria to explore chemically their environment.

Due to the small size of bacterial cells, the typical Reynolds number around a swimming bacterium is $Re\approx 10^{-4}$ and as a result, the fluid flow obeys the incompressible Stokes equations, which are time-reversible. The ability of bacteria to reorient requires them to break the time-reversibility constraint (i.e. the scallop theorem purcell1977 ), which is enabled by a short flexible hook ($\approx 60$ nm in length) that connects the rotary motor to the semi-rigid flagellar filament. The motor/hook/filament complex is known as a flagellum. The flexibility of the hook allows it to behave as a universal joint samatey2004 and is essential for flagellar bundle brown2012 . During the swimming motion, the flexible hook can buckle, providing a rich spectrum of swimming behaviours shum2012 ; son2013 ; nguyen2017 ; nguyen2018 ; jabbarzadeh2018 .

These mechanical structures together with the randomly distributed rotary motors of a peritrichous bacterium raise a fundamental question: In which direction does the cell move for a given motor configuration? The dynamics of a cell resulting from the rotating propulsion of multiple rotating objects (the flagellar filaments) could be referred to as $N$-flagella problem in reference to the classical $N$-vortex problem on the dynamics of point vortices newton2001 . The flagellar morphologies, including pitch size and the radius of the helical filament, fall onto one of a small number of polymorphic shapes, which have been characterised experimentally spagnolie2011 while the flagellar length can vary greatly within species and cell populations. Despite these variations, a coherent flagellar bundle is expected to be maintained generically by long-range hydrodynamic interactions watari2010 ; kanehl2014 .

Recently, Riley et al. riley2018 demonstrated theoretically that the swimming of peritrichous bacteria is enabled by an elastohydrodynamic instability. Modelling a bacterium as propelled by multiple rigid helical filaments spatially distributed around the cell body and connected to it by linear torque springs, they showed that the coupling between the flagellar propulsive force pushing the cell body and the hydrodynamic forces resulting from the swimming motion could result in an elastohydrodynamic instability of the hook and lead to bacterial swimming towards a preferred direction. This motility transition, demonstrated numerically, was also explained by a theoretical model of a bacterium propelled by two rod-like flagella connected to the spherical cell body from opposite sides riley2018 . The critical value predicted by the linear stability analysis was in good agreement with the full computational simulations, suggesting that the swimming direction of peritrichous bacteria might be set by the stability of an equilibrium distribution of flagellar filaments.

The aim of the current paper is to formalise this physical result mathematically and to derive rigorously the elastohydrodynamic motility transition theory for a cell with an arbitrary number of flagella of arbitrary shape using linear stability analysis. In the theoretical model of Ref. riley2018 , the rod-like flagella only generate forces, which the issue of torque generation was not considered. Since flagellar propulsion is generated by the rotation of helical filaments, a generation of torque is inevitable and we include it in this paper. The number of flagella is set to any integer $N$ greater than 1 and we will consider equilibria and its stability for the spatial distribution of the flagella.

The paper is organised as follows. Section II is devoted to the theoretical formulation of the problem. In Section IIA, we formulate the bacterial motility problem for a cell with $N$ flagella, considering both force and torque balance for the cell and each flagellum. In Section IIB, we proceed to simplify the problem by focusing on the case of identical flagellar filaments generating axisymmetric propulsion. The latter property holds for a helical filament assuming that the timescale of the bending and rotation of the flagella are well separated, which is verified in practice. The typical flagellar rotation frequencies are $\approx 100$ Hz and sufficiently faster than the typical frequency of cell rotation, $\approx 10$ Hz, that we may approximate the flagellar propulsion by its time-average. We then consider the case where the $N$ flagella are symmetrically distributed in a plane around the cell body and we focus on its linear stability from the equilibrium configuration, which is presented in In Section IIC. The following sections are devoted to discussions of the results of the linear stability analysis. We first neglect small chirality effects in order to simplify the system, and in Section III we start our analysis with the case of $N=2$ flagella, which is found to be different from the general case with $N>2$ due to the symmetry of the flagellar distribution. The general case is then discussed in Section IV, where we start with the examples of $N=3$ and $N=4$ before deriving rigorously the general stability results. Finally, we reincorporate the effects of chirality in Section V where solve the full problem.

II Mathematical model of peritrichous bacterium

II.1 Equations of motion

In this first section, we describe the force and torque balance equations, together with the torque balance of each elastic spring and hydrodynamic drag on each flagellar filament, in order to formulate a linear problem of $2N+6$ dimensions for bacterial motility propelled by rotating flagella.

We consider a swimming bacterium located in a Newtonian fluid of constant dynamic viscosity $\mu$. The cell body is assumed to be a sphere of radius $R$. We denote the centre of the sphere by $\bm{X}$ and its orientation by a unit vector $\bm{e}$. The bacterium is assumed to possess $N$ rigid flagellar filaments connected to the cell body at their base. The direction for the axis of each filament is measured by the unit vectors $\bm{e}^{(i)}$ for $i$-th flagellum ($i=1,2,\cdots,N$), as schematically illustrated in Fig. 1. We label each flagellar filament using the arc length, $s\in[0,L^{(i)}]$, measured from the flagellum-cell body connection (i.e. the location of the motor), where $L^{(i)}$ is the length of the $i$-th flagellum, whose shape is determined by its tangent vector $\bm{t}^{(i)}(s)$.

II.1.1 Force balance for a whole cell

We first consider the force balance equations for the entire cell, which, in the absence of inertia, state that the sum of the hydrodynamic forces on the cell body and all flagella must add up to zero.

Let $\bm{U}$ and $\bm{\Omega}$ be the linear and angular velocities of the cell, respectively. Neglecting hydrodynamic interactions between the cell body and the flagellar filaments, the hydrodynamic drag on the cell body, $\bm{F}_{body}$, is simply given from the Stokes law, by $\bm{F}_{body}=C_{D}\bm{U}$, where $C_{D}=-6\pi\mu R$ and hereafter we non-dimensionalise the length scale using $R=1$.

The hydrodynamic force on the $i$-th flagellum, $\bm{F}^{(i)}$, is obtained from the resistive-force theory of slender filaments, which predicts that the hydrodynamic force on a small segment of the flagellum, $d\bm{F}$, is linearly related to its local velocity relative to the background fluid, $\bm{u}^{(i)}$, as

[TABLE]

where $C_{t}$ and $C_{n}$ are the negative drag coefficient constants whose values depend on the flagellar slenderness parametergray1955 ; lauga2009 . Here, the local velocity $\bm{u}^{(i)}$ is the sum of the velocities of the cell body and the flagella. Introducing the angular velocity of each flagellum as $\bm{\omega}^{(i)}$ (Fig. 1), we have

[TABLE]

where

[TABLE]

is the position along the flagellar segment at arc length $s$. From (2), integrating over the entire flagellar filament leads to the hydrodynamic force on the $i$-th flagellum, with denoting the drag coefficient tensor in (1) by $\mathsf{C}_{i}$, given by

[TABLE]

Here, we have introduced skew-symmetric matrices $\mathsf{A}_{i}$ and $\tilde{\mathsf{A}}_{i}$ whose components are given, respectively, by

[TABLE]

where $\varepsilon_{pqr}$ is the Levi-Civita symbol and the Einstein summation convention is used over the repeated indices ($p,q,r=1,2,3$). The final force balance equations obtained by summing up the forces as

[TABLE]

II.1.2 Torque balance for a whole cell

Similarly to the above arguments, we now derive the expressions for the torque balance for an entire cell at the centre of the cell body. We again neglect hydrodynamic interactions between the cell and flagella and model hydrodynamic drag on each flagellar filament at the level of resistive force theory (RFT).

The torque acting on the spherical cell body is given by $\bm{M}_{body}=C_{R}\bm{\Omega}$, where the resistance coefficient is $C_{R}=-8\pi\mu R^{3}$. The hydrodynamic torque on a segment of a flagellum is given by $d\bm{M}^{(i)}=(\bm{x}^{(i)}+\bm{\xi}^{(i)})\times\,d\bm{F}^{(i)}$, which yields the total torque expression after an integration over the flagellum as

[TABLE]

where he superscript, $T$, indicates the transpose of a matrix. Note that in order to obtain that expression we neglected the torque arising from local rotation of the flagellar filaments around their centreline, which is typically orders of magnitude smaller than (7) in the experimental limit where the radius of the helical centreline is much larger than the thickness of the filament.

The overall torque balance equation for the whole cell is then written

[TABLE]

II.1.3 Torque balance for each flagellum

We proceed to consider the torque balance relation for each flagellum, which experiences both hydrodynamic torque and an elastic spring restoring torque at the flagellum-cell body junction ishimoto2017 ; riley2018 . The hydrodynamic torque follows from the discussion above and we now consider the torque balance at the flagellum-cell body junction point. The torque acting on a segment of a flagellum is then given by $d\tilde{\bm{M}}=\bm{\xi}^{(i)}\times d\bm{F}^{(i)}$ and integrating over the flagellum yields the hydrodynamic torque,

[TABLE]

As a model for the elastic hook, a linear spring torque is assumed to be present at the junction point. Let $\kappa^{(i)}>0$ be the spring constant and let the relative angle difference from the initial orientation be denoted by $\theta^{(i)}$. The elastic torque on each flagellum can be then written as

[TABLE]

where $\bm{e}^{(i)}_{\perp}$ is the unit vector perpendicular both to the initial and current flagellar orientation vectors. In the later part of this manuscript, we will assume that the flagellar orientations initially coincide with the outward normal $\bm{n}^{(i)}$ (Fig. 1), but the formulation here does not necessarily assume this initial condition.

Hereafter, we assume that each flagellum is rotated at a fixed rate and thus the torque balance equation for each flagellum is obtained as the instantaneous balance

[TABLE]

where we note that only the torque balance perpendicular to the vector $\bm{e}^{(i)}$ is considered, where $\mathsf{P}_{i}=(\mathsf{1}-\bm{e}^{(i)}\bm{e}^{(i)})$ is the projection onto the plane perpendicular to $\bm{e}^{(i)}$.

Note that if, alternatively, one was to model swimming as induced by motors rotating at fixed torque shum2012 ; kanehl2014 , we would need an additional term in the torque balance equation, which would take the form $\tilde{\bm{M}}^{(i)}+\bm{M}^{(i)}_{elast}+\bm{M}^{(i)}_{motor}=\bm{0}$, where now it is the moment $\bm{M}^{(i)}_{motor}$ which has a fixed value. In what follows, we focus on the rotation-given problem, aiming to generalise the results predicted by the previous theoretical model riley2018 . The theoretical extension to the torque-given problem with $N=2$ flagella is presented in Appendix B where we highlight the similarities and differences between the two models.

II.2 Identical and axisymmetric flagellar propulsion

From (4), (7) and (9), we obtain the governing equations in a matrix form,

[TABLE]

where $A$ is a square matrix of order $2N+6$, given by

[TABLE]

We next decompose the flagellar rotation velocity vector into the components due to the flagellar rotation and that due to the flagellar bending, $\bm{\omega}^{(i)}=\bm{\omega}^{(i)}_{t}+\bm{\omega}^{(i)}_{n}$, with $\bm{\omega}^{(i)}_{t}=\left(\bm{\omega}^{(i)}\cdot\bm{e}^{(i)}\right)\bm{e}^{(i)}$. In the rotation-given problem, $\omega^{(i)}_{t}=\bm{\omega}^{(i)}\cdot\bm{e}^{(i)}$ is a fixed value. For more convenience, we introduce $3\times 3$ matrices to simplify the linear equations (12) as

[TABLE]

In order to proceed, we assume that all $N$ flagella are identical and that they generate axisymmetric propulsion around their long axes. This assumption is obviously satisfied for axisymmetric flagellar shapes such as rods, but is also valid for asymmetric shapes, including helices, provided that the time-scale of flagellar rotation is much smaller than that of elastic bending (i.e., $|\bm{\omega}_{t}|\gg|\bm{\omega}_{n}|$). In practice, this assumption of axisymmetric propulsion is satisfied for swimming bacteria lauga2006 . Indeed, the typical flagellar rotation frequency of flagellar filaments for E. coli is $\approx 100$ Hz, which is faster than the typical frequency of cell rotation, $\approx 10$ Hz. We can therefore approximate the propulsion by a rotating helical flagellar filament by its time-averaged contribution, which is axisymmetric along the helix axis.

We introduce a reference frame, flagellum-fixed frames and the rotation matrix, mapping the laboratory reference frame, $\{\bm{e}_{x},\bm{e}_{y},\bm{e}_{z}\}$, onto the flagellum-fixed coordinates, $\{\bm{e}^{(i)}_{x},\bm{e}^{(i)}_{y},\bm{e}^{(i)}_{z}\}$ and denoted by $\mathsf{R}_{i}$, as schematically shown in Fig. 1. The axisymmetric conditions are satisfied if the quantity is unchanged under rotation around $\bm{e}^{(i)}_{z}$, where we choose the flagellum-fixed coordinates such that the $z$ axis ($\bm{e}^{(i)}_{z}$) coincides with $\bm{e}^{(i)}$ (See also Fig 1).

For any second-rank tensor written in the reference-frame as $\mathsf{K}$, its expression in the body-fixed frame, $\mathsf{K}_{i}$, is given by

[TABLE]

We next assume that the tensor is axisymmetric around the vector $\bm{e}^{(i)}=\bm{e}^{(i)}_{z}$. This can be expressed by requiring invariance under rotation by any angle $\alpha\in[0,2\pi)$ around the axis, i.e.

[TABLE]

where $\mathsf{R}\left(\alpha;\bm{e}^{(i)}\right)$ denotes the rotation matrix of angle $\alpha$ around $\bm{e}^{(i)}$. If the tensor $\mathsf{K}_{i}$ is represented as a matrix, the axisymmetric property leads to the matrix form,

[TABLE]

and the four upper-left entries can be decomposed into the symmetric and skew-symmetric parts,

[TABLE]

where $K_{s}$ and $K_{ss}$ are constants. From the matrix form in (17), note that we readily obtain the commutation relation,

[TABLE]

The tensors $\tilde{\mathsf{A}}_{i}$ and $\mathsf{C}_{i}$ depend only on the flagellar shape and if the flagellum shape is axisymmetric, these tensors are also axisymmetric. As a result the tensor $\mathsf{K}^{(i)}_{FF}=\int\tilde{\mathsf{A}}^{T}_{i}\cdot\mathsf{C}_{i}\cdot\tilde{\mathsf{A}}_{i}\,ds$ is found to be axisymmetric, since $\mathsf{K}^{(i)}_{FF}=\int\mathsf{R}\left(-\alpha;\bm{e}^{(i)}\right)\cdot\tilde{\mathsf{A}}^{T}_{i}\cdot\mathsf{C}_{i}\cdot\tilde{\mathsf{A}}_{i}\cdot\mathsf{R}\left(\alpha;\bm{e}^{(i)}\right)\,ds=\mathsf{R}\left(-\alpha;\bm{e}^{(i)}\right)\cdot\mathsf{K}_{FF}^{(i)}\cdot\mathsf{R}\left(\alpha;\bm{e}^{(i)}\right)$. If instead the flagellar filaments are not rigorously axisymmetric, notably if they are helical, then as long as the rotation velocity around their long axis is sufficiently large compared with the velocity scale during bending, we can approximately replace the second-rank tensors by their time-averages, and the same arguments thus follow.

Using the relation (19), equation (14) can be simplified to

[TABLE]

Similar arguments for $\mathsf{K}^{(i)}_{FF}$ enable us to show that the symmetric tensors $\mathsf{K}^{(i)}_{TF}=\mathsf{K}^{(i)T}_{FT}$ are also axisymmetric around the axis $\bm{e}^{(i)}$. However, the tensors $\mathsf{K}^{(i)}_{RF}=\mathsf{K}^{(i)T}_{FR}$ are not necessarily axisymmetric, since $\mathsf{A}_{i}$ are not axisymmetric except when each $\bm{e}^{(i)}$ is parallel to $\bm{x}^{(i)}$. We then decompose $\mathsf{A}_{i}=\mathsf{A}^{\prime}_{i}+\tilde{\mathsf{A}}_{i}$. Noting that $\mathsf{A}^{\prime}_{i}$ is independent of the flagellar shape, we rewrite as $\mathsf{K}_{RF}^{(i)}=\mathsf{A}_{i}^{{}^{\prime}T}\cdot\mathsf{K}_{TF}^{(i)}+\mathsf{K}_{FF}^{(i)}$ so that the tensors $\tilde{\mathsf{A}}_{i}$, $\mathsf{K}^{(i)}_{TF}$ and $\mathsf{K}^{(i)}_{FF}$ are axisymmetric around the vector $\bm{e}^{(i)}$. For convenience, we write $\mathsf{A}_{i}^{{}^{\prime}T}\cdot\mathsf{K}_{TF}^{(i)}=\mathsf{K}_{RF}^{{}^{\prime}(i)}$.

We next employ the assumption that all flagellar filaments have identical propulsion, which guarantees that the axisymmetric tensor, $\mathsf{K}$, can be expressed as $\mathsf{K}=\mathsf{R}_{i}\cdot\mathsf{K}^{(0)}\cdot\mathsf{R}_{i}^{-1}$, where the tensor $\mathsf{K}^{(0)}$ is common to all filaments. We also write the rotation of the flagellar angular velocity in the flagellum-fixed frame as $\bm{\omega}^{(i)}_{t}=\mathsf{R}_{i}\cdot\tilde{\bm{\omega}}^{(i)}_{t}$; here $\tilde{\bm{\omega}}^{(i)}_{t}$ is parallel to $\bm{e}_{z}$ and we can thus write $\tilde{\bm{\omega}}^{(i)}_{t}=\tilde{\omega}^{(i)}_{t}\bm{e}_{z}$, where the value of $\tilde{\omega}^{(i)}_{t}$ is prescribed. Similarly, the elastic torque in the flagellum-fixed frame is written as $\bm{M}^{(i)}_{elast}=\mathsf{R}_{i}\cdot\tilde{\bm{M}}^{(i)}_{elast}$. Since both $\mathsf{K}^{(i)}_{TF}$ and $\mathsf{K}^{(i)}_{FF}$ are axisymmetric, we can rewrite the linear equations (20) as

[TABLE]

where $\mathsf{P}_{0}$ is the projection onto the $x-y$ plane, defined as $\mathsf{P}_{i}=\mathsf{R}_{i}\cdot\mathsf{P}_{0}\cdot\mathsf{R}_{i}^{-1}$. With additional rotational matrix in the bottom $N$ rows, this can be simplified to

[TABLE]

Here, we have used the fact that $\mathsf{P}_{0}$ commutes with $\mathsf{K}^{(0)}_{FT}$ and $\tilde{\mathsf{K}}^{(0)}_{FR}$ and that $\mathsf{P}_{0}\cdot\tilde{\bm{\omega}}_{n}^{(i)}=\tilde{\bm{\omega}}_{n}^{(i)}$. The first and second entries of the right-hand side of (22) give the effective force and torque for the entire body, corresponding to propulsive force and torques,

[TABLE]

For a symmetric configuration such that $\left(\sum_{i=1}^{N}\tilde{\omega}^{(i)}_{t}\mathsf{R}_{i}\right)\cdot\bm{e}_{z}=\bm{0}$, only the first term of $\bm{M}_{prop}$ contributes to the motion of the cell. If we further assume that the reference flagellar orientation is normal to the sphere surface, i.e. $\bm{e}^{(i)}_{z}(t=0)=\bm{n}^{(i)}$, we readily obtain $\tilde{\bm{M}}^{(i)}_{elast}=\bm{0}$ and $\bm{M}_{prop}=\bm{0}$. Therefore, we have stationary solutions satisfying $\bm{U}=\bm{\Omega}=\tilde{\bm{\omega}}_{n}^{(i)}=\bm{0}$ for the symmetric-filament configurations.

From the axisymmetric propulsion of the flagellar filaments and using (17)-(18), we may obtain the general form of the tensors

[TABLE]

where the constants, $k_{C},K_{C},k_{T},k_{F}$, are all negative due to the negative definiteness of the Stokes resistance tensors and $K_{T}$ and $K_{F}$ are taken to be negative so that a positive value of $\tilde{\omega}_{t}$ generates flagellar force towards the cell body (as is the case for the majority of flagellated bacteria). From the definitions of $\mathsf{K}^{(0)}_{C}$ and $\mathsf{K}^{(0)}_{FF}$, these matrices are symmetric, thus the expressions in (25) follow. The expression for $\mathsf{K}_{TF}^{(0)}$ reflects the skew-symmetric property of $\tilde{A}_{i}$, and thus the diagonal component, $k_{D}$, is typically very small for swimming bacteria. The diagonal component of the matrix can be either positive or negative depending on the chirality of flagellar filaments. Using the Appendix in Ref. lauga2006 , $k_{D}$ can be estimated for a helical flagellum as $k_{D}\sim C_{N}bL\epsilon$ with $b$ is the diameter of the helix and $\epsilon=2\pi b/\lambda$, where $\lambda$ is the pitch of the helix. The ratio of diagonal to off-diagonal component is thus given by $k_{D}/k_{T}\sim\epsilon(b/L)\sim 10^{-2}$, always a small number for a bacterium such as E. coli lauga2006 . Note that $k_{D}$ is identically zero for rod-like flagella.

We now proceed with the linear stability problem neglecting the chirality effects (i.e. setting $k_{D}=0$) and will then incorporate chirality back in Section V.

II.3 Linear stability

In this section, we assume the in-plane flagellar configuration is that of a regular polygon, and then formulate the linear stability problem around the equilibrium. The example case of $N=6$ is schematically shown in Fig. 2(a) with notation. As shown above, with a symmetric flagellar configuration we obtain stationary solutions for the motion of a bacterium with $N(\geq 2)$ identical axisymmetric flagella. We assume that the angular rotations of each flagellum, $\bm{\omega}^{(i)}_{t}$, are identical and that each flagellar filament is connected perpendicularly to the cell body surface and located in one plane that we define as $x-y$ plane (see Fig. 2a). When the flagellar configuration forms an $N$-sided regular polygon, the bacterial motion is found to be stationary.

We now consider the linear stability around this equilibrium, where the orientation of $i$th flagellum is given by rotation of angle $\Theta_{i}$ from the $+x$ axis in the $x-y$ plane where $\Theta_{i}=2\pi(i-1)/N$ (see Fig. 2b). The disturbance from the equilibrium is denoted by two angles for each flagellum. We introduce the in-plane and out-of-plane angle displacements for the $i$th flagellum as $\theta_{i}$ and $\phi_{i}$, as shown in Fig. 2(a, b). In the linear stability regime, both angles are assumed to be sufficiently small, i.e. $|\theta_{i}|,|\phi_{i}|\ll 1$. We write the flagellar bending rotation as $\tilde{\bm{\omega}}^{(i)}_{n}=(d\theta_{i}/dt)\bm{e}^{(i)}_{x}+(d\phi_{i}/dt)\bm{e}^{(i)}_{y}$, assume a constant flagellar rotation around the axisymmetric axis, $\tilde{\bm{\omega}}^{(i)}_{t}=\omega_{0}$, and also assume that the torque spring constants are identical for all flagella, i.e. $\kappa^{(i)}=\kappa$.

We then proceed to obtain the expressions for the cell dynamics around the equilibrium configuration. The rotation matrix, $\mathsf{R}_{i}$ is obtained by combining the two rotations,

[TABLE]

where $\mathsf{R}(\phi_{i};\bm{e}_{y})$ denotes the rotation matrix with angle $\phi_{i}$ around the orientation $\bm{e}_{y}$. Using the expressions (25), we directly obtain the effective force,

[TABLE]

and similarly the effective torque,

[TABLE]

We then proceed to compute the force and torque generated by the in-plane flagellar bending, $\tilde{\bm{\omega}}_{n}^{(i)}=(d\theta_{i}/dt)\bm{e}_{x}$. The force generated by each flagellar filament is $\mathsf{R}_{i}\cdot\mathsf{K}_{TF}^{(0)}\cdot\tilde{\bm{\omega}}_{n}^{(i)}\simeq k_{T}(d\theta_{i}/dt)(-\sin\Theta_{i}\bm{e}_{x}+\cos\Theta_{i}\bm{e}_{y})$. Similarly, the torque generated can also be computed, using the decomposition $\mathsf{K}_{RF}^{(i)}\cdot\mathsf{R}_{i}\cdot\tilde{\bm{\omega}}_{n}^{(i)}=\left(\mathsf{A}_{i}^{{}^{\prime}T}\cdot\mathsf{R}_{i}\cdot\mathsf{K}_{TF}^{(0)}+\mathsf{R}_{i}\cdot\mathsf{K}_{FF}^{(0)}\right)\cdot\tilde{\bm{\omega}}_{n}^{(i)}=(k_{T}+k_{F})(d\theta_{i}/dt)\bm{e}_{z}$.

The force generated by the out-of-plane flagellar bending, $\tilde{\omega}_{n}^{(i)}=(d\phi_{i}/dt)\bm{e}_{y}$, is also computed as $\mathsf{R}_{i}\cdot\mathsf{K}_{TF}^{(0)}\cdot(d\phi_{i}/dt)\bm{e}_{y}=k_{T}(d\phi_{i}/dt)\mathsf{R}_{i}\cdot\bm{e}_{x}\simeq k_{T}(d\phi_{i}/dt)\bm{e}_{z}$ . Similarly, we obtain the torque generated by the bending, $(d\phi_{i}/dt)\bm{e}_{y}$, as $\mathsf{K}_{RF}^{(i)}\cdot\mathsf{R}_{i}\cdot(d\phi_{i}/dt)\bm{e}_{y}=(d\phi_{i}/dt)(\mathsf{A}_{i}^{{}^{\prime}T}\cdot\mathsf{R}_{i}\cdot\mathsf{K}_{TF}^{(0)}+\mathsf{R}_{i}\cdot\mathsf{K}_{FF}^{(0)})\cdot\bm{e}_{y}=(d\phi_{i}/dt)(k_{T}\mathsf{A}_{i}^{{}^{\prime}T}\cdot\mathsf{R}_{i}\cdot\bm{e}_{x}+k_{F}\mathsf{R}_{i}\cdot\bm{e}_{y})=(d\phi_{i}/dt)(k_{T}+k_{F})(\sin\Theta_{i}\bm{e}_{x}-\cos\Theta_{i}\bm{e}_{y})$.

A lengthy but straightforward calculation then leads to the matrix expression for the motion around the equilibrium, captured by a $2N+6$ dimensional linear problem,

[TABLE]

where $\bm{\Phi}=(U_{x},U_{y},U_{z},\Omega_{x},\Omega_{y},\Omega_{z},\dot{\theta}_{1},\cdots,\dot{\theta}_{N},\dot{\phi}_{1},\cdots,\dot{\phi}_{N})^{T}$ and

[TABLE]

The detailed derivations leading to the expressions for the square matrix $\mathcal{A}$ can be found in Appendix A and the results are

[TABLE]

where the constants in the matrix are

[TABLE]

and the $N$ dimensional vectors are given by

[TABLE]

and we used $\bm{1}_{N}$ and $\bm{0}_{N}$ to denote the identity and zero matrices of order $N$. In the following sections, we will consider the stability of this linear system.

III Instability of bacteria with $N=2$ flagella

We start with analysing the $N=2$ case. Substituting $N=2$ in the expressions (29)- (31), we obtain a 10-dimensional linear system (Fig. 3). With the geometric symmetry of the system, we readily find that the $x$ component of the effective force and torque vanish and thus $U_{x}=\Omega_{x}=0$ follow. The system is then reduced to an 8-dimensional linear system. Introducing the variables $\theta_{+}=\theta_{1}+\theta_{2}$, $\theta_{-}=\theta_{1}-\theta_{2}$, $\phi_{+}=\phi_{1}+\phi_{2}$ and $\phi_{-}=\phi_{1}-\phi_{2}$, we can rewrite the linear system (29) into 4 blocks of 2$\times$2 matrices as

[TABLE]

where we have defined $C_{D}^{\prime}=C_{D}+2k_{C}$ and $C_{R}^{\prime}=C_{R}+2k_{C}+4k_{T}+2k_{F}$ and where the dot symbol indicates the time derivative of the angle variable.

Solving each 2-dimensional problem yields ordinary differential equations with respect to the four angle variables as

[TABLE]

where the expressions of $A_{T},A_{R},A_{RT}$ are given by

[TABLE]

and where the two determinants, $\Delta_{D}=C_{D}^{\prime}k_{F}-2k_{F}^{2}$ and $\Delta_{R}=C_{R}^{\prime}k_{F}-2(k_{T}+k_{F})^{2}$, are positive as a result of the negative-definiteness of the Stokes resistance matrices.

We first consider the simple case where $K_{F}=0$ where the flagellar filaments produce propulsion but no torque. This assumption follows for a rod-like active filament as in Ref. riley2018 . Under this assumption, the matrix (45) is diagonal, since $A_{RT}=0$ and with eigenvalues $A_{T}$ and $A_{R}$. The angles $\theta_{-}$ and $\phi_{+}$ correspond to the translation modes in the $y$ and $z$ directions respectively (Fig. 3a). The critical flagellar angular velocity, $\omega_{0T}$, above which the translation mode becomes unstable is given by

[TABLE]

which is positive. When $\omega_{0}>0$, the flagellar filaments exert propulsive forces pushing on the cell body, whereas they pull on the organism when $\omega_{0}<0$. We therefore obtain that the translation instability only occurs for flagellar filaments in the pushing mode and with sufficiently large propulsive magnitude (or, for a fixed propulsion, with a sufficiently flexible hook). This can be compared to the arguments put forward in Ref. riley2018 , where the force by a rod-like flagellum corresponds to the product $\omega_{0}K_{T}$. Note that here the translation can occur towards an arbitrary direction in the $y-z$ plane.

The eigenvectors corresponding to the eigenvalue $A_{R}$ are linear combinations of the angles $\theta_{+}$ and $\phi_{-}$ and characterise a rotation mode around $z$ and $y$ axis, respectively (Fig. 3b). Although the rotation mode has not been examined in Ref. riley2018 (that study assumed mirror-image symmetry in their theoretical description), the rotation instability can occur above a second critical flagellar angular velocity, $\omega_{0R}$, given by

[TABLE]

Note that this critical value is positive and therefore the rotation instability also occurs only when the flagellar filaments in the pushing mode.

In the case where $K_{F}\neq 0$, each flagellar filament generates both torque and force. Again, we obtain the same eigenvalues of the matrix, $A_{R},A_{T}$, as in the problem of the flagella without torque generation. Moreover, the eigenvectors associated with the eigenvalue $A_{R}$ are the same as above (pure rotation move). As a difference, however, the eigenvectors for the eigenvalue $A_{T}$ are combinations of the four angles and the induced motion is found to be a translation in $y-z$ plane combined with a rotation around the translation direction. This unstable mode corresponds to the case of bacterial flagellar bundling via the elastohydrodynamic instability as observed experimentally and reproduced using numerical simulations in Ref. riley2018 .

III.1 Most unstable mode

When the angular velocity of each flagellar filament is positive, the two modes can become unstable. However, only the most unstable mode from the linear stability theory is likely to be observed in practice and we now consider which one of the two modes becomes unstable first.

Let us denote the length of the identical flagella by $L$. Using dimensional analysis allows to obtain order-of-magnitude estimates for the dependence with $L$ of the constants from the flagellar shape, namely $k_{C}=O(L)$, $k_{T}=O(L^{2})$ and $k_{F}=O(L^{3})$. Noting that $C_{D}=-6\pi$ and $C_{R}=-8\pi$, we can then estimate the size of the ratio $r$ between the two critical angular velocities, $r\equiv\omega_{0T}/\omega_{0R}$, as $r\sim C_{D}/C_{R}=3/4$ when $L\ll 1$ and $r=O(L^{-1})$ when $L\gg 1$. In both limits, we see that the critical value for the translation instability is smaller and thus it is the one which would be observed.

For further discussions in the intermediate region of $L$, we consider rod-like flagellar filaments of radius $d$ and length $L$. We introduce a positive constant $c$ such that $C_{n}=-4\pi\mu c$ is the normal drag coefficient, $c$ is given by $c=(\log(2L/d)+0.5)^{-1}$ lauga2009 . Within resistive force theory, the coefficients are given by $k_{C}=C_{n}L$, $k_{T}=(C_{n}/2)L^{2}$ and $k_{F}=(C_{n}/3)L^{3}$ and we can compute numerically the ratio $r$ for different flagellar parameters, $c$ and $L$, or $d$ and $L$, with results plotted in Fig. 4.

In Fig. 4(a), iso-values of the frequency ratio $r$ is first shown for different values $c$ and $L$ and we find that in the intermediate region $L\sim 1$ the rotation mode can be more unstable, although the translational mode is always more unstable above a critical value of $c\approx 0.4$. The same plot is then shown with the flagellar radius $d$ in the horizontal axis of Fig. 4(b), indicating that the rotation instability would occur first in a robust range of flagellar radius if $L\sim 1$.

Typical sizes of the cell body and the flagellar filaments of E. coli are $\approx 1~{}\mu$m and $\approx 10~{}\mu$m lauga2009 , leading to a large non-dimensional flagellar length, $L\approx 10$, and thus the translation mode is predicted to be the experimentally-observable one. Similarly, the typical size of the dimensionless flagellar radius and the RFT coefficient prefactor are given by $d\approx 0.02$ and $c\approx 0.13$, respectively, and a cell equipped with short flagellar filaments could therefore, in theory, undergo the rotational instability without net locomotion.

We next consider a helical flagellum whose shape is characterised by the helix angle $\Psi$ between the local flagellar tangent vector, $\bm{t}^{(i)}$ and the axis of the helix, $\bm{e}^{(i)}_{z}$. The tangent vector can be expressed by $\bm{t}^{(i)}=\cos\Psi\bm{e}^{(i)}_{z}+\sin\Psi(-\sin\beta\bm{e}^{(i)}_{x}+\cos\beta\bm{e}^{(i)}_{y})$ with the angle $\beta$ in the range $\beta\in[0,2\pi)$. For the flagellar filament moved along the $\bm{e}^{(i)}_{x}$ axis with velocity $u$, direct computations from (1) gives the local drag force,

[TABLE]

where $\gamma=C_{t}/C_{n}$ is the ratio of the tangential and normal coefficients from resistive force theory. If the timescale of flagellar rotation is sufficiently faster than that of bending, we can approximate the local force by its time-averaged value, which is obtained averaging over the angle parameter $\beta$. We thus obtain $d\bm{F}\cdot\bm{e}^{(i)}_{x}=C^{\prime}_{n}u$, with the effective normal drag coefficient $C^{\prime}_{n}=-4\pi\mu c^{\prime}$ as

[TABLE]

enabling us to follow the same arguments for the rod-like flagella by simply replacing $C_{n}$ by $C^{\prime}_{n}$, or $c$ by $c^{\prime}$. Typical values of $\gamma$ and $\Psi$ for E. coli are $\gamma\approx 0.7$chattopadhyay2006 and $\Psi\approx 30$ deg.spagnolie2011 , resulting in the effective drag coefficient, and thus the effective value of $c$, to be similar to the case of a rod ($c^{\prime}\approx 0.13$). As before, the translation mode is therefore more unstable for a bacterium with typical sizes, though there is a region where the rotation mode can be more unstable when $L\sim 1$. Note however that this assumes that the helical structure of the filament is neglected, a modelling assumption which is corrected in Section V.

IV Instability of bacteria with $N\geq 3$ flagella

We now proceed to study the linear stability problem in the general case of $N$ flagellar filaments. When $N\geq 3$, the drag coefficients (32) - (37) can be simplified to

[TABLE]

where we have used the equalities,

[TABLE]

which are satisfied when $N\geq 3$. We first address the $N=3$ and $N=4$ cases to allow for intuition on the mathematical structure of the solution, before proceeding to the general $N$ case.

IV.1 Cell with $N=3$ flagella

When $N=3$, the number of the angular variables is 6, which is equal to the number of degrees of freedom for the rigid motion of the whole cell. As in the previous section, we partially diagonalise the matrix (31) by introducing the angle variables, $\tilde{\theta}_{1}=\theta_{2}-\theta_{3}$, $\tilde{\theta}_{2}=2\theta_{1}-\theta_{2}-\theta_{3}$, and $\tilde{\theta}_{3}=\theta_{1}+\theta_{2}+\theta_{3}$ for the in-plane angles and as $\tilde{\phi}_{1}=\phi_{2}-\phi_{3}$, $\tilde{\phi}_{2}=2\phi_{1}-\phi_{2}-\phi_{3}$, and $\tilde{\phi}_{3}=\phi_{1}+\phi_{2}+\phi_{3}$ for the out-of-plane angles.

We then obtain a linear system decomposed into 6 two-by-two block matrices in which one degree of freedom for rigid motion is paired with an angle variable. For the translation and rotation in $x$ direction, the matrices are given by

[TABLE]

from which we obtain the ordinary differential equations for the linear stability in the form,

[TABLE]

where $A_{T}=\Delta_{T}^{-1}((3/2)|K_{T}k_{T}|\omega_{0}-|C_{D1}|\kappa)$, $A_{R3}=\Delta_{R3}^{-1}((3/2)|K_{T}(k_{T}+k_{F})|\omega_{0}-|C_{R1}|\kappa)$ and $A_{TR3}=(3/2)\Delta_{R3}^{-1}|K_{F}(k_{T}+k_{F})|\omega_{0}$, with the determinants $\Delta_{T}=C_{D1}k_{F}-(3/2)k_{T}^{2}$ and $\Delta_{R3}=C_{R1}k_{F}-(3/2)(k_{T}+k_{F})^{2}$. From (58), we can read off directly the eigenvalues for the linear stability with the rigid motion in $x$ direction as $A_{T}$ and $A_{R3}$. The eigenvectors are pure rotation around the $x$ axis ($A_{R3}$ eigenvalue) and translation along the $x$ axis accompanied by rotation around the $x$ axis ($A_{T}$).

In a similar manner, we can derive the eigenvalues for linear stability associated with the rigid motion in the $y$ direction. Noting that $C_{D1}=C_{D2}$ and $C_{R1}=C_{R2}$, we have the similar equation to (58),

[TABLE]

which yields the same eigenvalues $A_{T}$ and $A_{R3}$ as in the linear stability in $x$ direction.

For the translation and rotation in the $z$ direction, the matrices are now given by

[TABLE]

yielding the equation for the linear stability,

[TABLE]

where $A_{T3}=\Delta_{T3}^{-1}(3|K_{T}k_{T}|\omega_{0}-|C_{D3}|\kappa)$, $A_{R3}=\Delta_{R}^{-1}(3|K_{T}(k_{T}+k_{F})|\omega_{0}-|C_{R1}|\kappa)$ and $A_{TR3}=3\Delta_{R}^{-1}|K_{F}(k_{T}+k_{F})|\omega_{0}$, with the determinants $\Delta_{T3}=C_{D3}k_{F}-3k_{T}^{2}$ and $\Delta_{R}=C_{R3}k_{F}-3(k_{T}+k_{F})^{2}$. The system in (62) provides the eigenvalues of the linear stability associated with the rigid motion in $z$ direction as $A_{R}$ and $A_{T3}$, with corresponding eigenvectors of pure rotation around the $z$ axis and combined translation and rotation along the $z$ axis, respectively.

In summary, we obtain 6 eigenvalues, $A_{T}$, $A_{T}$, $A_{T3}$, $A_{R}$, $A_{R3}$, $A_{R3}$, each of which is accompanied by a rigid-motion mode for the whole cell.

IV.2 Cell with $N=4$ flagella

In the case of a cell equipped with $N=4$ flagella, the number of angle variables exceeds the numbers of degrees of freedom for a rigid-body motion of the cell body. We again introduce new angle variables to decompose the square matrix of order $2N+6$ into smaller size systems as $\tilde{\theta}_{1}=\theta_{2}-\theta_{4}$, $\tilde{\theta}_{2}=\theta_{1}-\theta_{3}$, $\tilde{\theta}_{3}=\theta_{1}+\theta_{2}+\theta_{3}$ and $\tilde{\theta}_{4}=\theta_{1}-\theta_{2}+\theta_{3}-\theta_{4}$, with similar combinations for the out-of-plane angle variables. From these changes of variables, we obtain the 6 matrices associated with the rigid motion and the angle variables, $\tilde{\theta}_{i}$ and $\tilde{\phi}_{i}$ for $i=1,2,3$. In turn, we obtain the same form of the linear ordinary differential equations for the linear stability as of (58), (59) and (62), for the $x$, $y$ and $z$ directions, respectively. However, the values of the matrix entries are now given by $A_{T}=\Delta_{T}^{-1}(2|K_{T}k_{T}|\omega_{0}-|C_{D1}|\kappa)$, $A_{R3}=\Delta_{R3}^{-1}(2|K_{T}(k_{T}+k_{F})|\omega_{0}-|C_{R1}|\kappa)$ and $A_{TR3}=(2\Delta_{R3}^{-1}|K_{F}(k_{T}+k_{F})|\omega_{0}$, with the determinants $\Delta_{T}=C_{D1}k_{F}-2k_{T}^{2}$ and $\Delta_{R3}=C_{R1}k_{F}-2(k_{T}+k_{F})^{2}$ for the instabilities in $x$ and $y$ directions and $A_{T3}=\Delta_{T3}^{-1}(4|K_{T}k_{T}|\omega_{0}-|C_{D3}|\kappa)$, $A_{R3}=\Delta_{R}^{-1}(4|K_{T}(k_{T}+k_{F})|\omega_{0}-|C_{R1}|\kappa)$ and $A_{TR3}=4\Delta_{R}^{-1}|K_{F}(k_{T}+k_{F})|\omega_{0}$, with the determinants $\Delta_{T3}=C_{D3}k_{F}-4k_{T}^{2}$ and $\Delta_{R}=C_{R3}k_{F}-4(k_{T}+k_{F})^{2}$ for the instability in the $z$ direction.

The remaining two angular degrees of freedom are diagonalised as

[TABLE]

which leads the two negative eigenvalues for the linear stability, $\kappa/k_{F}<0$. This, in turn, indicates that the instability can occur only if accompanied by the rigid-body motion of the cell. Here we note that the angle variables, $\tilde{\theta}_{4}$ and $\tilde{\phi}_{4}$, do not generate any forces and torque. Inspecting the definition of $\tilde{\theta}_{4}$, we see that the force from the angles $\theta_{1}-\theta_{2}$ is canceled by that generated by $\theta_{3}-\theta_{4}$ and the torque from the part $\theta_{1}+\theta_{3}$ cancels that from $\theta_{2}+\theta_{4}$.

IV.3 General N case

The two simple examples above enable us to now characterise the instabilities in the general $N$ case. Specifically, we expect that the instabilities are associated with the rigid-body translation and rotation of the cell body even when $N\geq 5$. We introduce the new angle variables as found in the expressions of $\bm{F}_{prop}$ (27) and $\bm{M}_{prop}$ (28)

[TABLE]

We first consider the mode associated with the translation and rotation along the $x$ axis. The resulting $2\times 2$ matrices are

[TABLE]

from which we obtain two eigenvalues,

[TABLE]

where we have introduced the determinants $\Delta_{T}=C_{D1}k_{F}-(N/2)k_{T}^{2}$ and $\Delta_{R3}=C_{R1}k_{F}-(N/2)(k_{T}+k_{F})^{2}$ which are positive due to the negative-definiteness of the resistance matrices. The expressions (68) and (69) are similar to the results obtained for $N=3$ and $N=4$.

As expected by symmetry, the modes associated with the translation and rotation along $y$ are similar and the eigenvalues are the same as in the $x$ direction, i.e. $A_{T}$ and $A_{R3}$.

We then proceed to investigating the modes along the $z$ axis and obtain the $2\times 2$ matrices

[TABLE]

This system has eigenvalues

[TABLE]

with determinants $\Delta_{T}=C_{D1}k_{F}-Nk_{T}^{2}$ and $\Delta_{R3}=C_{R1}k_{F}-N(k_{T}+k_{F})^{2}$ The eigenvalues, (72) and (73), again reproduce the results of $N=3$ and $N=4$.

The eigenvectors obtained so far are associated with three pure rotational modes and three combined translation/rotation modes along the same axis. The remaining degrees of freedom associated with the other $2N-6$ angular variables do not affect the stability of the cell, which can be summarised to the following statement:

Theorem IV.1

The linear system (12) includes the six modes associated with a rigid-body motion, with eigenvalues $A_{T}$, $A_{T}$, $A_{T3}$, $A_{R}$, $A_{R3}$, $A_{R3}$. The remaining $2N-6$ degrees of freedom all generate identical negative eigenvalues, $\kappa/k_{F}<0$.

In order to complete the proof of the above statement, we need to rearrange the angle variables so as to diagonalise the remaining $2N-6$ degrees of freedom. We first prepare linearly independent $N$ in-plane angle variables, $\tilde{\theta}_{i}$, for $i=1,\cdots,N$, in which the angle variables $\tilde{\theta}_{1},\tilde{\theta}_{2},\tilde{\theta}_{3}$ defined in (64) are included. When $N\geq 4$, the diagonalisation can then be achieved if we pick the remaining angle variables so that they do not generate any effective force and torque. This is possible for an arbitrary $\tilde{\theta}_{i}$ for $i\geq 4$, since we can add $\tilde{\theta}_{1},\tilde{\theta}_{2},\tilde{\theta}_{3}$ without disobeying the linear independence property for the set of the angle variables. This argument can also be applied to the out-of-plane variables and the remaining $2N-6$ eigenvalues are found to be all identical and negative, $\kappa/k_{F}<0$.

IV.4 Most unstable mode

We obtain $2N$ eigenvalues for the linear system, six of which can be positive. In this section, we examine the nature of the most unstable mode, which is the one expected to be relevant in an experiment. From equations (68), (69), (72) and (73), we obtain critical angular velocities, above which the system becomes linearly unstable, as

[TABLE]

These values are all positive, indicating that the system is stable in the case where the flagellar filaments pull on the cell body and that the instability can occur only when flagella push. From equations (53) and (54), we readily obtain the comparison between the critical values for the in-plane and out-of-plane instabilities as

[TABLE]

which indicate that the most unstable mode is either the translation towards the $z$ axis (with combined rotation around the same direction) or pure rotation around the $z$ axis.

Notably, the expressions of the critical values include those of the $N=2$ case, (49) and (50). However, due to the symmetry of the system, the eigenvalues are degenerated and the relation (76) becomes equalities $\omega_{0T3}=\omega_{0T}$ and $\omega_{0R}=\omega_{0R3}$. Nonetheless, the expressions of $\omega_{0T3}$ and $\omega_{0R}$ in (74)-(75) can be obtained when we substitute the constants $k_{C}\mapsto(N/2)k_{C}$, $k_{T}\mapsto(N/2)k_{T}$ and $k_{F}\mapsto(N/2)k_{F}$ in (49)-(50). Thus, the rod-like flagellar model examined for $N=2$ can also be applied in the general $N$ case, and replacing $c\mapsto(N/2)c$ we obtain the same plots for the ratio of the two critical values, $r$, as in Fig. 4(a). This indicates that the increase in the value of $N$ can remove the possibility of $r>1$, and therefore the translation instability would always be the most unstable mode for pushing flagella.

We note that simulation results with $N=4$, which corresponds to a typical number of flagella for E. coli berg2003 , showed translation in one direction along with rotation around the same axis riley2018 , a result consistent with our stability analysis. We also note that the critical values (74) and (75) do not depend on the value of the torque generated by the flagellar rotation, $K_{F}$, but they only depend on the force generated by each rotating flagellum. In contrast, the value of $K_{F}$ appears in the eigenvector of the translation mode and cells with large $K_{F}$ undergo rapid rotation.

V Small chirality effects

The analysis in the previous sections neglected the diagonal components of the matrix $\mathsf{K}^{(0)}_{TF}$. In this section, we reincorporate the diagonal terms $k_{D}$ and solve the resulting linear stability problem. Since the diagonal components are small compared to the off-diagonal terms, we may treat the full problem as a perturbation from the results obtained in the previous sections.

As in the previous sections, the linear stability problem is written $\mathcal{A}\bm{\Phi}=\bm{b}$ where straightforward calculations now lead the matrix (see also Appendix A)

[TABLE]

with new components arising from the chirality given by

[TABLE]

and new $N$-dimensional vectors

[TABLE]

With the use of the new angle variables (64) - (65), we may again reduce the linear problem into $4\times 4$ blocks associated with the $j$th ($j=1,2,3$) component of the force and torque, while the remaining degrees of freedom yield only negative eigenvalues for the stability problem. The explicit form along the $x$ direction can be computed as

[TABLE]

Inverting the matrix on the left hand side in the previous equation leads to the following linear ordinary differential equations

[TABLE]

for which we can characterise the linear stability as in the previous sections. Using the small parameter $\delta$ to measure the relative magnitude of the chirality effects in the resistance matrix, $\delta\sim k_{D}/k_{F}\sim 10^{-2}$ for the typical parameters of E. coli bacteria, we may expand the matrix in Eq. (85) at first order in $\delta$ as $\tilde{A}_{T}=A_{T}+\delta A^{\prime}_{T}$, $\tilde{A}_{R3T}=\delta A^{\prime}_{R3T}$, $\tilde{A}_{TR3}=A_{TR3}+\delta A^{\prime}_{TR3}$, $\tilde{A}_{R3}=A_{R3}+\delta A^{\prime}_{R3}$. This results in similar eigenvalues and eigenvectors for the system with changes of magnitude of order $\delta$. However, as a result of the coupling term $\tilde{A}_{R3T}$, the pure rotation modes no longer exist and all eigenmodes are now associated with translation in $x$ direction (note that since the changes due to chirality are small, former rotation modes are still dominated by rotation and generate small net locomotion).

The analysis along the $y$ and $z$ directions are similar and we obtain perturbed eigenvectors and eigenvectors from those obtained in the previous sections. Except for the 6 angle variables, all other modes continue to have stable eigenvalues $\kappa/k_{F}(<0)$ as shown in the previous section, which can be summarised into the following statement:

Theorem V.1

The linear stability problem characterised by matrix (77) includes six modes associated with rigid-body motion, and the remaining $2N-6$ degrees of freedom all generate identical negative eigenvalues, $\kappa/k_{F}<0$.

The elastohydrodynamic instability can occur if the flagellar filaments push the cell body, accompanied by the rigid motion of the whole cell, with critical angular velocities perturbed from those obtained in the previous section, $\omega_{0T},\omega_{0T3},\omega_{0R},\omega_{0R3}$, since $\delta$ is small. All unstable modes now lead to transition in one direction accompanied by rotation along the same direction. However, there are two types of modes for each direction, translation-dominated and rotation-dominated modes, the latter of which becomes pure-rotation modes when chirality effects are neglected.

VI Discussion

In this paper, we investigated theoretically the elastohydrodynamic stability problem of a model bacterium with multiple flagellar filaments rotated with prescribed frequencies. We assumed that the cell was equipped with $N$ identical flagella connected to a spherical body surface by a flexible elastic torque spring and that the flagella are initially arranged in a plane with equal angle intervals as to form a regular $N$-polygon. We first formulated the equations of motions of this system and showed that this configuration provides an equilibrium state where the cell body does not move.

We then proceeded to consider the linear stability problem in the case of negligible chirality in the flagellar filaments (active rods). When $N=2$, two modes are obtained (translation and rotation) which can be unstable when the flagella push on the cell body provided the magnitude of this pushing force exceeds a critical value (or, for a fixed propulsion magnitude, provided the hook is sufficiently flexible). The translation mode is the more unstable for the typical parameters of real bacteria and corresponds to the translation in one direction with rotating around the same axis. However, when the flagellar lengths are of the same order as the cell radius, $L\sim R$, the most unstable mode could be switched to the second mode where the cells rotate in place in the plane of the initial flagellar configuration with no associated translation.

We then extended our results to the general case of $N$ flagellar filaments, and we found that there are always only 6 modes which can be unstable, all of which are associated with rigid-body motion of the cell. The most unstable mode induces translation towards the direction perpendicular to the plane in which the $N$ flagella are initially arranged and is accompanied by rotation around the same axis. This analytical result is in agreement with numerical simulation with $N=4$ helical flagella riley2018 .

We finally reincorporated the chirality of the flagellar filaments which had been neglected in the previous sections. Chirality leads to small perturbations of the eigenvalues and eigenvectors for the linear stability and there are still only 6 possible unstable modes for cells with pusher flagella associated with rigid motion of the whole cell. The rotation modes are now accompanied by a small translation of the cell due to chirality-induced coupling.

The theoretical results in this paper and the presence of rotation-dominated modes thus imply that multi-flagellated peritrichous bacteria with shorter flagella could fail to swim efficiently. In contrast, for cells with typical flagellar length $L\sim 10$, the flagella produce a sufficient amount of propulsive thrust to lead to an instability to translation. A similar analysis could also be applied to synthetic particles propelled by bacterial flagella darnton2004 ; dileonardo2010 and to the dynamics of an ovum pushed by multiple spermatozoa ishimoto2017 . Note that the linear stability analysis performed in this paper can obviously not fully predict the nonlinear dynamics after the initial stages of the instability. Furthermore, as shown in Appendix B, if one considers instead the case of flagella rotated by a constant torque applied in the direction normal to the cell surface, all the possible unstable modes are accompanied by cell translation and pure-rotation modes disappear. This is in contrast to the case where the constant torque is applied along the long axis of the flagellar filament, for which the stability analysis coincides with the fixed-rotation case (see Appendix B). These results emphasise the complexity of the multi-flagellated swimming dynamics.

Using typical parameter values for E. coli ($R=1\mu$m, $L=10\mu$m, $N=4$) and the value of the viscosity for water ($\mu=10^{-3}$Pa$\cdot$s), we can estimate the critical flagellar rotation rate provided by our theory. The strength of the torque spring for an E. coli hook has been estimated to be in the range $\kappa\approx 2.9-8.7\times 10^{-21}$Nm using measured values fo the hook bending stiffness and length son2013 ; riley2018 , and we use the value $\kappa=5\times 10^{-21}$Nm for the following discussions. From the flagellar propulsion force used in Riley et al. riley2018 we have $K_{T}\approx 7.0\times 10^{-16}$N$\cdot$s. The drag coefficients are estimated as in Sec. III.1 with $c^{\prime}\approx 0.13$ and can be used to obtain an estimation of the critical rotation frequencies for a E. coli cell with $N=4$ flagellar filaments. The critical values for the translation modes are predicted to be $\nu_{0T}\approx 0.36$Hz and $\nu_{0T3}\approx 0.29$Hz, which are small compared with those for the pure rotation mode, $\nu_{0R}\approx 1.33$Hz and $\nu_{0R3}\approx 1.34$Hz. Since the flagellar filaments of real cells rotate much faster ($\nu\approx 100$Hz), the elastohydroynamic instability obtained in this paper is likely to be relevant to the locomotion of bacteria.

The model in our paper could be readily extended to the case of a spheroidal cell body if the case where the flagella are all initially arranged in the equatorial plane of the spheroid, and one would simply need to change the values of the drag coefficient $C_{D}$. When the cell body takes the shape of a prolate spheroid such as of E. coli, $C_{D1}/C_{D3}>2$ still holds and we obtain the same relation as (76). Other straightforward extensions include the situation in which the cell body is located near a planar infinite wall, a situation relevant to a sperm-egg cluster that tends to rotate without translation ishimoto2017 . The predominance of rotation could be rationalised in that case using lubrication theory lauga2006 , which shows that drag coefficients $C_{D1}$ and $C_{D3}$ diverge as the spherical cell body approaches the wall, while the value of the rotation drag coefficient $C_{R3}$ remains very close to that in the bulk. As a result, the rotation mode would become in that case more unstable than the translation mode. Further theoretical work would be required to extend to more general situations such as for example a non-spherical cell body, non-symmetric flagellar configurations, or non-identical flagella, emphasising the rich diversity of the $N$-flagella problem.

Data Accessibility. This paper has no additional data.

Authors’ Contributions. KI and EL designed the research and developed the mathematical model. KI analysed the model. KI and EL wrote the paper.

Competing Interests. The authors declare that they have no competing interests.

Funding. This project has received funding from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation programme (grant agreement 682754 to EL). KI is supported by MEXT Leading Initiative for Excellent Young Researchers (LEADER), JSPS KAKENHI (Grand Number JP18K13456) and JSPS Overseas Research Fellowship (29-0146).

Acknowledgements. We thank anonymous referees for helpful comments.

Appendix A. Derivations of matrices (31) and (77)

In this Appendix, we provide detailed derivations of the matrices (31) and (77) for the motion of the bacterium with $N$ flagella arranged in a regular polygonal manner.

The $N$ identical in-plane flagella with orientations $\bm{e}^{(i)}=\mathsf{R}(\Theta_{i};\bm{e}_{z})\cdot\bm{e}_{x}$ provide a stationary configuration. Considering small disturbances around the equilibrium, with angles $|\theta_{i}|,|\phi_{i}|\ll 1$, we directly compute the matrix entries noting that one only needs the leading order contributions for the linear stability.

We first consider $\mathsf{K}_{C}^{(i)}$, which is computed from (25),

[TABLE]

where the symbol $\simeq$ is used here to mean the leading order contribution. The expression for $K_{TT}$ is given by summation of this matrix over the indices $i$.

For the off-diagonal part, we need $\sum_{i=1}^{N}(\mathsf{K}_{C}^{(i)}\cdot\mathsf{A}_{i}^{\prime}+\mathsf{K}^{(i)}_{TF})$. From expression (86), we have

[TABLE]

and using (25) we obtain

[TABLE]

The summations in (87) and (88) give the expression for $\mathsf{K}_{TR}$ and its transpose just follows for $\mathsf{K}_{RT}$. The contributions of the $K_{T}$ terms can be neglected following the approximation showing that $k_{D}$ is negligible. For a helical filament, $K_{F}$ scales as $K_{F}\sim C_{N}b^{2}L$, and comparing it with the leading-order term we have the relative magnitude as $K_{F}/k_{F}\sim(b/L)^{2}\sim 10^{-3}$ using typical numbers for E. coli cells. Thus in the matrix (31) the $K_{F}$ term can be neglected if we neglect $k_{D}$.

For the expression of $\mathsf{K}_{RR}$, we need to calculate

[TABLE]

which is obtained by straightforward calculations as

[TABLE]

[TABLE]

[TABLE]

Using the equalities,

[TABLE]

and summing over the index $i$ completes the computations for the matrix entries.

Appendix B. Torque-driven motility of $N=2$ flagella

In this Appendix, we briefly consider the bacterial model with $N=2$ flagella in case where the torque, instead of the rotation, is prescribed for each flagellar filament. This will allow us to highlight the difference of the dynamics from the rotation-given problem. As in the main text, we assume identical flagella and axisymmetric propulsion. Using the same matrix form as (14), the torque-driven motility dynamics can be expressed as

[TABLE]

We use the decomposition of the flagellar rotation velocity vector, $\bm{\omega}^{(i)}=\bm{\omega}^{(i)}_{t}+\bm{\omega}^{(i)}_{n}$, and the commutative relations (19), to obtain the same form of the force and torque balance equations as (20), namely

[TABLE]

The torque balance equation for each flagellum is

[TABLE]

which can be rewritten, using the variables in the flagellum-fixed frame, as

[TABLE]

We next introduce the projection on to $z$ axis as $\mathsf{Q}_{0}=\mathsf{1}-\mathsf{P}_{0}$, and apply $\mathsf{P}_{0}$ and $\mathsf{Q}_{0}$ from the left side of the equation (97). From the projection onto the $x-y$ plane, we obtain a similar torque balance equation as in the bottom row of the equation (20), namely

[TABLE]

noting that the last term of the right-hand side is the only correction from the rotation-given problem.

The projection using $\mathsf{Q}_{0}$ provides the equations for the tangential flagellar rotation velocity,

[TABLE]

and we proceed to calculate the detailed expressions for the linear stability analysis around the equilibrium configuration with $N=2$.

We need to assume the exact form of the function $\tilde{\bm{M}}^{(i)}_{motor}$, and here we consider two different possibles for the constant torque: (i) constant torque applied along the flagellar orientation ($\bm{M}^{(i)}_{motor}=M_{0}\bm{e}^{(i)}$), and (ii) constant torque applied along the normal to the cell surface ($\bm{M}^{(i)}_{motor}=M_{0}\bm{n}^{(i)}$).

When $|\theta_{i}|,|\phi_{i}|\ll 1$, the leading-order value of the right-hand side of (99) is given by $-M_{0}\bm{e}_{z}$ in both torque models. The third term on the left-hand side of (99) is simply $\mathsf{K}^{(0)}_{FF}\cdot\tilde{\bm{\omega}}^{(i)}_{t}=K_{F}\omega_{0}\bm{e}_{z}$, and thus we can neglect $O(|\theta_{i}|,|\phi_{i}|)$ contribution in the first two terms of (99) in order to determine the leading-order term of $\omega_{0}$. With calculations similar to those in Appendix A, we obtain the $O(1)$ contribution from the first term as $\mathsf{Q}_{0}\cdot\mathsf{K}^{(0)}_{FT}\cdot\mathsf{R}^{-1}_{i}\cdot\bm{U}\simeq K_{F}(\cos\Theta_{i}U_{x}+\sin\Theta_{i}U_{y})\bm{e}_{z},$ which is however zero as a consequence of the fact that $U_{x}=\sin\Theta_{i}=0$ for the linear stability problem with $N=2$. The second term is calculated as $\mathsf{Q}_{0}\cdot\mathsf{R}^{-1}_{i}\cdot\mathsf{K}^{(i)}_{FR}\cdot\bm{\Omega}\simeq K_{F}(\Omega_{x}\cos\Theta_{i}+\Omega_{y}\sin\Theta_{i})\bm{e}_{z}$ and this is again found to be zero since $\Omega_{x}$ and $\sin\Theta_{i}$ are zero for the linear stability problem with $N=2$.

In summary, we obtain the expression for the flagellar rotation rate, $\omega_{0}$,

[TABLE]

and $\omega_{0}$ becomes positive when $M_{0}>0$. Equations (95) and (98) provide therefore a set of equations similar to (20) for the rotation-given problem. We note the presence of the additional term $\mathsf{P}_{0}\cdot\tilde{\bm{M}}^{(i)}_{motor}$ in equation (98). This correction term, however, vanishes for the torque model (i), and thus this torque-driven motility problem coincidences with the rotation-given motility when $N=2$.

We then consider the torque model (ii) where the constant torque is applied along the normal to the cell surface. In that case, the correction term contributes an external bending torque, since $\mathsf{P}_{0}\cdot\tilde{\bm{M}}^{(i)}_{motor}\simeq M_{0}(\theta_{i}\bm{e}_{y}-\phi_{i}\bm{e}_{x})$. Proceeding with the linear stability analysis, as in equations (41)- (44) we obtain 4 blocks of $2\times 2$ matrices, using the angle variables $\theta_{+}=\theta_{1}+\theta_{2}$, $\theta_{-}=\theta_{1}-\theta_{2}$ , $\phi_{+}=\phi_{1}+\phi_{2}$ and $\phi_{-}=\phi_{1}-\phi_{2}$,

[TABLE]

and

[TABLE]

Solving equations (101)-(102) with respect to $\theta_{+}$ and $\phi_{+}$, we obtain the linear ordinal differential equations,

[TABLE]

The diagonal components are

[TABLE]

are the same as in the rotation-given problem (46)- (47). In contrast, the off-diagonal components includes the corrections

[TABLE]

and therefore the stability characteristics could be different from the rotation-given problem.

When $M_{0}<0$, the flagellar filaments pull on the cell-body into fluid and we find that the eigenvalues are all negative since $A_{T},A_{RT},A_{R}<0$ and $A_{TR}>0$. Hence the dynamics is always linearly stable. In contrast, when $M_{0}>0$, the eigenvalues are still both negative for small values of $M_{0}$, but become positive when $M_{0}$ becomes sufficiently large, and there are two critical values of $M_{0}$ above which instability can occur. Both unstable modes, however, combine translation with rotation and the pure rotation modes disappear (in contrast to the rotation-given motility problem). The critical flagellar rotation velocities $\omega_{0}^{\ast}$ in this case lie between the critical values of the rotation-given problem (49)- (50) :

[TABLE]

Note that, by symmetry, the same eigenvalues follow for the linear ordinal differential equations obeyed by the variables $\theta_{-}$ and $\phi_{i}$ derived from equations (103)-(104).