Article Information

• PDF (1752 kb)
• Full Text with Large Figures

PubMed

• Articles by Charles W. Wolgemuth

Related Articles

• Equilibrium Shapes of Erythrocytes in Rouleau Formation

  Equilibrium Shapes of Erythrocytes in Rouleau Formation Biophysical Journal, Volume 84, Issue 3, 1 March 2003, Pages 1486-1492 Jure Derganc, Bojan BoŽič, Saša Svetina and Boštjan Žekš Abstract A well known physiological property of erythrocytes is that they can aggregate and form a rouleau. We present a theoretical analysis of erythrocyte shapes in a long rouleau composed of cells with identical sizes. The study is based on the area difference elasticity model of lipid membranes, and takes into consideration the adhesion of curved axisymmetric membranes. The analysis predicts that the erythrocytes in the rouleau can have either a discoid or a cup-like shape. These shapes are analogous to the discoid and stomatocyte shapes of free erythrocytes. The transitions between the discoid and cup-like shapes in the rouleau are characterized. The occurrence of these transitions depends on three model parameters: the cell relative volume, the preferred difference between the areas of the membrane bilayer leaflets, and the strength of the adhesion between the membranes. The cup-like shapes are favored at small relative volumes and small preferred area differences, and the discoid shapes are favored at large values of these parameters. Increased adhesion strength enlarges the contact area between the cells, flattens the cells, and consequently promotes the discoid shapes. Abstract | Full Text | PDF (153 kb)

• Bacterial Flagellar Microhydrodynamics: Laminar Flow over Co...

  Bacterial Flagellar Microhydrodynamics: Laminar Flow over Complex Flagellar Filaments, Analog Archimedean Screws and Cylinders, and Its Perturbations Biophysical Journal, Volume 85, Issue 3, 1 September 2003, Pages 1345-1357 Shlomo Trachtenberg, Dalia Fishelov and Matania Ben-Artzi Abstract The flagellar filament, the bacterial organelle of motility, is the smallest rotary propeller known. It consists of 1), a basal body (part of which is the proton driven rotary motor), 2), a hook (universal joint—allowing for off-axial transmission of rotary motion), and 3), a filament (propeller—a long, rigid, supercoiled helical assembly allowing for the conversion of rotary motion into linear thrust). Helically perturbed (so-called “complex”) filaments have a coarse surface composed of deep grooves and ridges following the three-start helical lines. These surface structures, reminiscent of a turbine or Archimedean screw, originate from symmetry reduction along the six-start helical lines due to dimerization of the flagellin monomers from which the filament self assembles. Using high-resolution electron microscopy and helical image reconstruction methods, we calculated three-dimensional density maps of the complex filament of H13-3 and determined its surface pattern and boundaries. The helical symmetry of the filament allows viewing it as a stack of identical slices spaced axially and rotated by constant increments. Here we use the closed outlines of these slices to explore, in two dimensions, the hydrodynamic effect of the turbine-like boundaries of the flagellar filament. In particular, we try to determine if, and under what conditions, transitions from laminar to turbulent flow (or perturbations of the laminar flow) may occur on or near the surface of the bacterial propeller. To address these questions, we apply the boundary element method in a manner allowing the handling of convoluted boundaries. We tested the method on several simple, well-characterized cylindrical structures before applying it to real, highly convoluted biological surfaces and to simplified mechanical analogs. Our results indicate that under extreme structural and functional conditions, and at low Reynolds numbers, a deviation from laminar flow might occur on the flagellar surface. These transitions, and the conditions enabling them, may affect flagellar polymorphism and the formation and dispersion of flagellar bundles—factors important in the chemotactic response. Abstract | Full Text | PDF (371 kb)

• Depolymerization-Driven Flow in Nematode Spermatozoa Relates...

  Depolymerization-Driven Flow in Nematode Spermatozoa Relates Crawling Speed to Size and Shape Biophysical Journal, Volume 94, Issue 10, 15 May 2008, Pages 3810-3823 Mark Zajac, Brian Dacanay, William A. Mohler and Charles W. Wolgemuth Abstract Cell crawling is an inherently physical process that includes protrusion of the leading edge, adhesion to the substrate, and advance of the trailing cell body. Research into advance of the cell body has focused on actomyosin contraction, with cytoskeletal disassembly regarded as incidental, rather than causative; however, extracts from nematode spermatozoa, which use Major Sperm Protein rather than actin, provide at least one example where cytoskeletal disassembly apparently generates force in the absence of molecular motors. To test whether depolymerization can explain force production during nematode sperm crawling, we constructed a mathematical model that simultaneously describes the dynamics of both the cytoskeleton and the cytosol. We also performed corresponding experiments using motile spermatozoa. Our experiments reveal that crawling speed is an increasing function of both cell size and anterior-posterior elongation. The quantitative, depolymerization-driven model robustly predicts that cell speed should increase with cell size and yields a cytoskeletal disassembly rate that is consistent with previous measurements. Notably, the model requires anisotropic elasticity, with the cell being stiffer along the direction of motion, to accurately reproduce the dependence of speed on elongation. Our simulations also predict that speed should increase with cytoskeletal anisotropy and disassembly rate. Abstract | Full Text | PDF (963 kb)

• …more

Copyright © 2008 The Biophysical Society. All rights reserved.
Biophysical Journal, Volume 95, Issue 4, 1564-1574, 15 August 2008

doi:10.1529/biophysj.107.118257

Biophysical Theory and Modeling

Collective Swimming and the Dynamics of Bacterial Turbulence

Charles W. Wolgemuth^,  

University of Connecticut Health Center, Department of Cell Biology and Center for Cell Analysis and Modeling, Farmington, Connecticut 06030-3505

Address reprint requests to Charles W. Wolgemuth, Dept. of Cell Biology, University of Connecticut Health Center, 263 Farmington Ave., Farmington, CT 06030-3505. Tel.: (860) 679-1655.

The general dogma of fluid dynamics defines two distinct dynamical regimes that can be distinguished by the dimensionless Reynolds number (Re), a ratio of inertia to viscosity. In terms of characteristic length and velocity scales, L and U, respectively, Re=UL/η, where η is the kinematic viscosity of the fluid. When Re<1, which occurs at low velocities and small length scales, gradients in the hydrostatic pressure are balanced by viscous shear, the flow is smooth and laminar, and vortices are seldom seen. At large velocities and large length scales (Re≫1), the viscosity is negligible compared to inertia. Therefore, gradients in the pressure act to accelerate the fluid. In this regime, sharp transitions in velocity are possible, which can be accompanied by vortices.

Bacterial swimming occurs in the low Reynolds number regime. For example, Bacillus subtilis is a rod-shaped cell that is roughly 4μm long and 0.7μm in diameter 1. At low concentrations, these cells swim at ∼20μm/s 2. The Reynolds number for a single cell is therefore ∼10⁻⁴. As the cell concentration increases, fluid flow between neighboring cells will feel drag due to the velocity difference between the cells and the fluid. Therefore, one expects that as the bacterial concentration increases, the effective Reynolds number for the system should decrease. However, experiments on the swimming behavior and fluid dynamics in regions of high concentrations (in excess of 10⁹ cells per cm³) of B. subtilis reveal very unexpected results 1,3. In these dense regions of bacteria, transient jet-like fluid motion occurs in which the speed of the jets is larger than the speed of an individual bacterium (∼50μm/s) 1. Furthermore, these jet-like regions are flanked by vortices that persist for timescales of ∼1s 1. The fluid dynamics on length scales of 100–200μm is reminiscent of a von Karman vortex street, which arises in fluids when the Reynolds number is >50 4. Due to the qualitative similarity to high Reynolds number fluid dynamics, these swimming patterns have been termed bacterial “turbulence”.

This fluid behavior seems counterintuitive. If the bacteria move independently of the fluid, then the fluid flow must interpenetrate the dense suspension of bacteria. Increasing the concentration of the bacteria should increase the drag between fluid and the bacteria, thereby increasing the effective viscosity of the solution and lowering the Reynolds number. Yet, at face value, the opposite effect is seen; increasing the concentration of bacteria leads to fluid flows that seem to have a much higher Reynolds number. This leads to a number of interesting physical questions. What physics is required to describe the dynamical behavior of these systems? Can an effective Reynolds number for this system be defined, and, if so, how? What sets the length and timescales for the jets and vortices? Why are the velocities that are observed larger than the swimming velocity of an individual bacterium?

Recent theoretical work sheds some light on these issues. In Simha and Ramaswamy 5, the fluid stress from a collection of force dipoles was added to the Navier-Stokes equations. Linear analysis of this model in the low Reynolds number regime showed that it was unstable under large wavelength disturbances and that fluctuation in the number density of particles is anomalously large. Other authors have also used nematic hydrodynamic equations coupled to local stress generation to derive continuum equations for active filament solutions 6,7. These articles have shown that contractile activity can cause dramatic differences in rheological behavior and can lead to modulated flowing phases and macroscopic phase separation. Simulations of discrete numbers of force dipoles 8 and extended dipole rods 9 have shown that vortices and jets can arise in these systems; however, these simulations neglected stochastic effects such as diffusion and run-and-tumble swimming of the bacteria.

In this work, I develop a two-phase model for the collective swimming of dense colonies of bacteria. This model treats the fluid and bacteria as independent, interpenetrating continuum phases. The bacteria are rod-like and therefore entropic effects will favor local alignment of the bacteria. Since the thrust on the bacteria comes from the flagellar bundle yet the drag on the bacteria acts on both the flagella and the cell body, the location of thrust is displaced from the center of drag, and the bacteria exert a dipole force on the fluid. Analysis of the run-and-tumble nature of bacterial swimming leads to a term that governs how tumbling can reorient the bacterial alignment. In this model, stochastic effects are treated using a mean-field approximation. For reasonable values of the model parameters, this model accurately reproduces the complex behavior observed in these systems and leads to a qualitative understanding of the length, time, and velocity scales that are present in the experiments.

The length scales over which vortical patterns and jets are observed in dense communities of swimming bacteria are much larger than the size of an individual bacterium. Therefore, on the scale of these flow patterns, it is reasonable to assume that the bacteria and the fluid can each be described as a continuum field. In this section, I develop a two-phase model to describe the average behavior of this system on these length scales. I begin by considering the swimming behavior of a single bacterium. Then, by defining the average density of the bacteria at a given location using a volume fraction, I construct energy and Rayleighan functionals to describe the dissipative dynamics of the two-phase system. The rod shape of the bacteria defines a local average orientation, and torque balance is used to define the dynamics of this bacterial orientation.

B. subtilis is a peritrichously flagellated bacterium with a roughly cylindrical cell body of radius, R∼0.4μm, and length, L∼4μm 1. Each bacterium has ∼10 flagella 10, which are helical filaments of length L_f∼10μm. The radius of a flagellum, R_f, is only 10nm, though, so the volume occupied by the flagella can be considered to be negligible compared to the cell body. Therefore, the volume of a single bacterium is V_b∼πR²L.

As a bacterium is an elongated object, its long axis defines an orientation. I define a unit vector, d, that lies along this direction. When all the flagella rotate in the same fashion (i.e., counterclockwise), the flagella bundle together. The coordinated rotation of the flagella exerts a propulsive force on the fluid. An equal but opposite force is exerted on the bacterium from the fluid. I assume that this thrust force lies along d with magnitude F₀ (Fig. 1). This force is the stall force, i.e., the force that would be required to stop a bacterium swimming in a stationary fluid.

Display large version of this figure

Figure 1

Schematic of a rod-shaped swimming bacterium in an external pressure gradient. The bacterium has length L and radius R. The bacterium exerts a force on the fluid F₀=2πη_fV₀d, where d is the direction of the long axis of the bacterium and V₀ is the free swimming velocity of the bacterium. The thrust of the bacterium is produced by rotation of the flagella, and the center of thrust is defined to be at x_CT. The fluid exerts a drag force, F_drag, on the bacterium. Because the bacterial body experiences drag but does not exert thrust, the center of drag, x_D, is not located at x_T but instead is positioned at x_D=x_T+bd.

As the Reynolds number is ≪1, inertial terms are negligible. Therefore, in a quiescent fluid, the drag force from the fluid must balance the thrust force so that the sum of the forces on the bacterium is zero. According to slenderbody hydrodynamics, the drag force on an elongated object is proportional to the velocity difference between the fluid and the swimming bacterium. Treating both the cell body and the flagella as filamentary objects and ignoring the pitch of the flagella, the drag coefficient for motion along the long axis of the bacterium is ζ_sh=2πη_fL/(ln(L/2R)+c)+2πη_fL_f/(ln(L_f/2R_f)+c) 11, where η_f is the fluid viscosity and c is a constant of order 1. Here the first term is the drag on the cell body and the second term is the drag on the flagella; I have treated the flagellar bundle as a single filamentary object. If an isolated bacterium swims with an average velocity V₀ in a stationary fluid, then F₀=ζ_shV₀.

Thrust force is localized at the flagella, whereas drag force is distributed over the entire bacterium. Therefore, the centers of thrust and drag are not coincident. I denote the centers of thrust and drag as x_T and x_D, respectively, and assume that these locations are separated by a distance b along the direction d: x_T=x_D−bd (Fig. 1).

In a fluid with an externally generated pressure gradient (i.e., a pressure gradient not generated by the swimming bacterium), the bacterium also experiences a force from the pressure field. Because the Stokes equations are linear, the pressure field can be broken into two components. One describes the pressure field generated by the motion of the bacterium, the force from which is included in the slenderbody drag coefficient. The other piece is the external pressure field. In a constant external pressure gradient, the total force on the bacterium from this component of the pressure field is

(1)

(2)

Other forces can be added to this description to take into account Brownian motion or other random processes as well as interactions between other bacteria. However, for my purposes here, I would like to describe a dense colony of bacteria for which it will be easier to use a continuum description.

It is convenient to describe a continuum field of bacteria using the volume fraction, ϕ_b, which is the fraction of space filled by bacteria in a given volume, V. The fluid volume fraction is ϕ_f=1−ϕ_b. I choose a simplified method for describing the volume fraction. Treating the center of drag of the ith bacterium, x_D,i, as a point source for volume, the volume fraction at point x can be defined as

(3)

Using this definition of ϕ_b, it is straightforward to calculate the average of bacterial variables inside V. For example, I define the average velocity of the bacteria at point x as

(4)

(5)

(6)

(7)

Viscous fluid dynamics is inherently dissipative. In the low Reynolds number regime, dissipative forces (and forces that arise from energy production) are balanced by conservative forces. An elegant method for handling dissipative dynamics is to define an energy functional, E, that describes the conservative forces and a Rayleigh dissipative functional, R, to describe the dissipative forces. This technique has been used previously in the context of two-phase flows to describe polymer solutions 12. An alternative method for deriving two-phase flow equations directly averages the dynamic equations and leads to an equivalent formulation 13. The two-phase formalism has been used extensively. For example, it has been used to describe the dynamics of the eukaryotic cytoskeleton 14,15, the swelling of porous media 16, tissue deformation 17, and the nonlinear swelling of polymer gels 18.

For the collective swimming of bacteria, I define two components of the energy. The first is due to random motions of the bacteria and fluid, which leads to entropic mixing. I use the Flory-Huggins mixing energy 19,

(8)

Since the bacteria experience thermal fluctuations and may also execute run-and-tumble swimming, I use the fluctuation-dissipation theorem to define an energy scale, ζ_shD_b, for the bacterial mixing, where D_b is a diffusion coefficient that takes into account thermal and swimming stochasticity. k_BT is thermal energy, and V_f is the size of a water molecule. The second term in this equation is the fluid mixing entropy. For ϕ_b≪1, a redefinition of the mechanical pressure can account for this term. As I will be considering situations when ϕ_b<0.2, I will neglect it in the rest of my calculations.

The second component of the energy comes from the constraint that in a given volume the sum of the bacterial and fluid volume fractions must sum to 1. Defining a Lagrange multiplier, P, the constraint energy is

(9)

There are two sources of dissipation. First, the fluid is viscous. In addition, I also allow for the possibility that bacterial collisions and interactions between the flagella may lead to an effective bacterial viscosity, η_b. The rod shape of the bacteria can also give rise to a number of other dissipative terms that are not considered here (see Sonnet et al. 20). The Rayleigh dissipative functional for these viscous interactions is

(10)

Viscous drag between the fluid and bacteria is the second form of dissipation. I use the fact that the drag force is proportional to the velocity difference between the bacteria and fluid and define a drag coefficient, ζ. Therefore, the drag dissipation functional is

(11)

To estimate the value for ζ, I assume that the total drag per volume is equal to the sum of the drags on each bacterium individually; i.e., I neglect the effect of hydrodynamic interactions between the bacteria on ζ. Therefore, ζ=ϕ_bζ_sh/V_b∼2η_fϕ_b(L+L_f)/R²L, where I have left off the logarithmic dependence in the denominator. This value is close to the first order approximation derived from the Darcy permeability for flow between a collection of spheres 21. For simplicity, I have chosen a scalar, ζ, which assumes that drag on the bacteria does not depend on whether the relative velocity of the bacteria is along d or perpendicular to it; however, slenderbody hydrodynamics predicts that the drag coefficient for motion perpendicular to the long axis is twice that for motions along the long axis.

The rotation of the flagella inputs energy into the system. Therefore, I also add an energy production term to the dissipative functional. The total force produced at a given position x comes from all the bacteria whose center of thrust lies at point x. Since it is possible for bacteria whose centers of drag are located at different locations to have the same center of thrust, the total force at a given location is a sum over all bacteria whose centers of thrust coincide with the point of interest; i.e., the total force produced at x is

(12)

(13)

The dynamic equations for the bacterial and fluid phases, respectively, are derived from these functionals using a variational principle:

(14)

(15)

(16)

In a similar fashion, I can derive the force balance equation for the fluid phase:

(17)

(18)

Using Eq. (18) and integrating the thrust force by parts, Eq. (17) simplifies to

(19)

The total force on the two-phase system is the sum of Eqs. (16):

(20)

(21)

To close the dynamical equations, I need to describe the dynamics for the average bacterial orientation. Since the bacteria are rod shaped, at high concentrations it is expected that there will be an entropically driven tendency for the bacteria to align, much like what occurs in nematic liquid crystals 22,23. Nematic liquid crystal theory defines a critical concentration, ϕ_cr∼R/L_eff, above which this alignment will occur. Here, L_eff takes into account the length of the bacterium as well as any effects due to the flagella. Therefore, I expect that ϕ_cr is ∼1%. The concentrations at which these turbulent flow patterns occur have been achieved only with swimming bacteria. Therefore, whether or not nonmotile bacteria at these concentrations have a tendency to align is not known.

In two dimensions, the direction d can be defined by the angle θ measured from the x direction: d=cosθx+sinθy. The entropic tendency to align gives rise to an effective elasticity, such that there is an energetic cost for splay, twist, and bend variations in the direction d. The magnitude of the elastic restoring force is set by the three Franck elastic constants K₁ (splay), K₂ (twist), and K₃ (bend). In two dimensions, the twist energy is zero. I choose K₁=K₃=K. If the bacteria behave as rigid rods, K=cɛϕ²L²/R³24,25,26,27, where c is a constant of proportionality and ɛ is an energy scale. For rotational diffusion driven by thermal energy, ɛ=k_BT, with k_B being Boltzmann's constant and T the temperature. In addition, this same elasticity can give rise to additional Ericksen stresses in the fluid equation 22,28; however, I assume that these terms are negligible in comparison to the dipole stresses induced by self-propulsion, which should be true at long wavelengths in the ordered state.

If the rotational diffusion is driven by the run-and-tumble behavior of the bacteria, then the fluctuation-dissipation theorem suggests that ɛ should be ∼ζ_shD_b. For this work, I choose the constant of proportionality to be in the range (0.2–0.4)/π²26. Drag due to rotations of the bacteria is balanced by the nematic, elastic torque. Rotational velocity of the fluid, ω_f=(1/2)∇×v_f, can exert torque on the bacteria. I assume a dimensionless drag coefficient for this torque, ζ_f. For real situations, I expect that all the bacteria are not perfectly aligned along the direction d. Therefore, there is a spread in the angular distribution of the bacteria. Run-and-tumble behavior of the bacteria will enhance this angular distribution. In the Appendix A Changes in angle due to run-and-tumble diffusion,Appendix B Numerical implementation of the model equations, I show that a spread in the angular distribution leads to diffusional rotation of the bacteria. A simplified form of the full hydrodynamics of nematics (see, for example, DeGennes and Prost 22 and Yue et al. 29) is used to define the rotation of the orientational field, d:

(22)

(23)

The closed system of equations that define the two-phase collective swimming model are Eqs. (6).

The model Eqs. (6) form a complex system of partial differential equations. To begin, it is useful to define characteristic scales for the problem to explore the behavior that is expected from this system. Experiments on these systems 1,3 have measured the characteristic length, time, and velocity scales. In Dombrowski et al. 1, the characteristic time, τ_c, and length, L_c, were defined using the temporal and spatial correlations of the velocity field, respectively. These two scales correspond roughly to the persistence time and size of the vortices. In addition, the maximum velocity in the jets was used to define the characteristic velocity, V_c. It was observed that V_c∼L_c/τ_c, as is expected.

In the model presented here, orientation of the bacteria is maintained by the nematic elasticity, the strength of which is set by K. A parameter with dimensions of Length²/Time (i.e., a diffusion-like coefficient) can be defined as K/ζ_r. Using this parameter, I can define a Peclet number as the dimensionless ratio of velocity to diffusion, Pc∼ζ_rV_cL_c/K. The size of the vortices should occur at around the length where advection begins to dominate the nematic elasticity, Pc>1. Therefore, I expect the characteristic length to be L_c∼K/ζ_rV_c.

The dipole stress from the bacteria that is exerted on the fluid is primarily balanced by the fluid viscosity. Assuming that the bacterial volume fraction is ϕ₀ and balancing the viscous force with the dipole force, , or V_c∼ϕ₀V₀bL_c/R². Combining these estimates leads to the relations

(24)

I solved the model equations on a two-dimensional, 400μm×400μm, periodic domain using COMSOL Multiphysics (COMSOL, Burlington, MA). Details of the numerical method are given in the Appendix A Changes in angle due to run-and-tumble diffusion,Appendix B Numerical implementation of the model equations. For my simulations, I used the parameter values that are given in Table 1 and varied V₀, K/ζ_r, R, and ζ_f. All simulations used an initial condition with ϕ(x,t)=0.15 and θ(x,t)=0 plus a small random perturbation.

When V₀=5μm/s, K/ζ_r=10⁻⁵cm²/s, R=0.7μm, and ζ_f=1, the uniform initial condition quickly (on the order of seconds) becomes unstable and develops into a velocity field filled with jets and vortices accompanied by a nonuniform volume fraction. This turbulent behavior is initiated by the volume fraction developing into a roll-like pattern (Fig. 2). Soon after this, the velocity field begins to form vortices (Fig. 2). Within a few seconds, a persistent, turbulent-looking dynamical system is observed (Fig. 3).

Display large version of this figure

Figure 2

Onset of bacterial turbulence. A time series showing the bacterial volume fraction (color map) and the fluid velocity field (arrows). The simulations begin with a uniform distribution of bacteria with d directed primarily along the x direction with a small random perturbation applied (t=0s). After 3s, the density begins to vary into a skewed roll-like pattern. After 6s, the rolls have broken apart and vortices have formed. By 9s, a “turbulent”-looking velocity field has emerged with high velocity jets flanked by vortices. The density fluctuations have also begun to coarsen. The range of the color map is different from frame to frame to emphasize the density variations. The parameters that were used are given in Table 1 with V₀=5m/s, K/ζ_r=10⁻⁵cm²/s, R=0.7μm, and ζ_f=1.

Display large version of this figure

Figure 3

Time series showing the evolution of the fluid velocity field and bacterial density variations at times well beyond the onset of the turbulent behavior. The color map shows the bacterial volume fraction, and the arrows show the velocity field. The same parameters were used as in Fig. 2.

In agreement with the experimental measurements 1,3, the vortices and jets in this system are transient (Fig. 3) with a timescale on the order of 1 s. To make a quantitative comparison to the experiments, I calculate the temporal correlation, and the spatial correlation, where I is averaged over orientations of r. I find correlations similar to those observed in the experiments 1 (Fig. 4). By fitting J(t) to an exponentially decaying function, I find a timescale of 0.54s. The anticorrelation length scale from I(|r|) is estimated by eye to be ∼70μm. In addition, I find velocities as large as 100μm/s, which is significantly larger than the free swimming velocity that I have used for this simulation.

Display large version of this figure

Figure 4

Temporal (a) and spatial (b) correlation functions of the fluid velocity field for V₀=5μm/s, K/ζ_r=10⁻⁵cm²/s, R=0.7μm, and ζ_f=1.

For this simulation, I observe a volume fraction variation of ∼35%. The magnitude of these variations have not been measured, but density fluctuations are observed. Another feature observed in these simulations is that the velocity typically points perpendicular to the gradients of the volume fraction.

To test the qualitative understanding that was presented above, a number of simulations were run that varied systematically the parameters V₀, K/ζ_rϕ, R, and ζ_f. The scaling relations 24 suggest that τ_c should be inversely proportional to V₀, independent of K/ζ_rϕ, and vary as R². In Figure 5ace, τ_c is plotted with respect to these functions of the parameters. The scaling of τ_c with V₀ and R is seen to obey the scaling relations quite well over the biologically realistic range that was considered. From the simulations, τ_c was found to scale with K^1/2, which was not accounted for in the scaling relations.

Display large version of this figure

Figure 5

Dependence of the characteristic time and velocity on V₀, K/ζ_rϕ, and R. All plots are done to show the correspondence to the scaling relations given in Eq. (7). (a, c, and e) τ_c versus 1/V₀, (K/ζ_rϕ)^1/2, and R². (b, d, and f) V_c versus , (K/ζ_rϕ)^1/2, and 1/R. The solid lines correspond to the best fit straight line through the simulation points. All other parameters are fixed and given in Table 1.

The characteristic velocity is expected to be proportional to and to (K/ζ_rϕ)^1/2 and be inversely proportional to R. In Figure 5bdf, V_c is plotted with respect to these functions of the parameters. For all three parameters, the scaling relations for V_c are observed to hold. At small values of R, though, I did notice some deviation from the expected scaling relations. I could not find a quantitative method for extracting the characteristic length from the spatial correlation of the fluid velocity field. However, I assume that L_c is proportional to V_cτ_c, which is consistent with the derivation of the scaling relations.

I also ran simulations to test the effect of the fluid torque on the bacteria, the magnitude of which is determined by the parameter ζ_f. As the bacterial force and dipole stress act along the direction vector of the bacteria, the fluid velocity and bacterial velocity are strongly entrained to point in the same direction, even in the absence of this torque. The fluid torque acts to accentuate this orientational entrainment. Therefore, I assumed that the dynamics would not be strongly affected by the value of this parameter. Varying ζ_f between 0 and 1 produced a relative change in τ_c of ∼40% and a relative change in V_c of ∼20%. I also tested the dependence of the pattern formation on the magnitude of the angular spread in the orientational dynamics. I found that decreasing the value of c in Eq. (23) produced a relative increase in τ_c of ∼50% and a relative decrease in V_c of ∼10%. Therefore, I do not expect run-and-tumble dynamics to strongly influence the pattern formation and fluid dynamics in these systems.

In this work, I have developed a two-phase model for the swimming of dense colonies of bacteria. Using physiological parameters, the model recreates accurately the flow fields that are observed in communities of B. subtilis. The model is necessarily complex, taking into account swimming and alignment of the bacteria, the dipole stress that the bacteria exert on the fluid, torque from the fluid, and run-and-tumble diffusion of the bacteria. The derived and validated scaling arguments, though, shed light on the salient parameters that govern the behavior of the macroscopic flow field. In essence, the flow velocities are set by a balance between the viscosity of the fluid and the dipole stress from the bacteria, and the vortex length scale is set by a competition between advection of the bacteria and nematic alignment. The timescale for the transient nature of the flow is given by the ratio of the vortex size to the velocity.

The turbulent dynamics that are observed in this system are the long time outcome of the instability predicted in Simha and Ramaswamy 5 for suspensions of self-propelled particles with orientational order. Recent numerical studies of suspensions of discrete, self-propelled, dipole rods have also observed this instability as well as the subsequent turbulent motions 9. The model presented here advances this previous work by doing a more general derivation of the equations of motion, performing the force balance on the fluid and bacteria separately. The numerical solution of these equations shows that for physiologically realistic parameters, the model accurately describes the experimental data. In addition, the scaling arguments that are derived and validated here provide an intuitive view of the turbulent behavior in this system.

This model produces testable predictions for the dependence of the flow fields on the parameters of the system. Primary features that the characteristic scales are sensitive to are the dimensions of the individual bacteria. The length and radius of the bacterium and the length of the flagella set the drag, the nematic elasticity, and the size of the dipole stress. Therefore, similar experiments done with different species of bacteria can provide a method for testing the validity of this model. In addition, the scaling relations suggest how the characteristic scales should depend on the average density of the bacteria. Systematically varying the density of the bacteria would also provide a test for this model.

It is interesting that the fluid flows that arise in this system are so similar to high Reynolds number flows. Specifically, the fact that von Karman street patterns with jets flanked by offset vortices are observed in flows with a Reynolds number >50 leads one to wonder what physical parameter in this low Reynolds number regime may act like a large Reynolds number. One possible solution based on this model is the effective Peclet number that was derived for the bacterial alignment, Pc∼ζ_rV_cL_c/K=η_fL²ϕ₀V_cL_c/6R²K. Using V_c=50μm/s, L_c=50μm, and the parameters given in Table 1, Pc∼200.

Beyond the intrinsic interest in exploring an intriguing and somewhat counterintuitive biophysical phenomenon, understanding the behavior of dense colonies of bacteria has industrial and biomedical applications. The development of biofilms, which can cause infection and corrode equipment, are dependent upon the aggregation of bacteria. Also, in dead-end filtration, where water with suspended bacteria is forced through a microfilter, the bacterial concentration in the filter can reach 10⁶–10⁸cm⁻³30. How bacterial swimming affects this process has not been looked at yet.

Appendix A

Changes in angle due to run-and-tumble diffusion

During bacterial swimming, reversal of the flagellar motor causes the bacteria to reorient 31. These reversals happen on average every 1s 31. Due to these random reorientations, I expect that the bacteria in a unit volume will not all have exactly the same direction, d, but rather the bacteria will have an angular distribution about this preferred direction. To analyze how run-and-tumble behavior can affect the dynamics of the preferred direction, I define the number density, n(x,ψ,t), which is the number of bacteria at point x at time t aligned with angle ψ relative to the x direction. The average bacterial orientation at position x, θ, is, therefore,

(A1)

The number density of bacteria will obey the continuity equation,

(A2)

where v=v(x,ψ,t). Taking the time derivative of Eq. (A1) and using Eq. (A2), I find that

(A3)

where I have defined the total number density of bacteria at To continue, I break the bacterial velocity into two pieces, v_b=V+u, V is independent of ψ, and u depends on ψ and is assumed to be directed along the bacterial direction, . Plugging this into Eq. (A3) leads to

(A4)

To make further progress in analyzing Eq. (A4), I assume that n is normally distributed about θ with a standard deviation of σ, . Using the identity

(A5)

and defining the average velocity

(A6)

Eq. (A4) becomes

(A7)

Inspection of Eq. (A7) reveals that the average alignment of the bacteria is advected by the total average bacterial velocity, V+〈u〉, and also reorients due to the spread in orientations from tumbling. To gain insight into the behavior of the reorientational term, I approximate the Gaussian function as a step function that equals 1 for θ−σ<ψ<θ+σ and is 0 everywhere else, and I choose u=V₀, the free swimming velocity of the bacteria. This assumption gives that N=2n₀σ and the dynamic equation for θ is

(A8)

where c=σsinσ. For Escherichia coli, a single tumbling event reorients the bacterium by 58° 32, giving c=0.86. In this work, I choose c=1.

It is interesting to look at the behavior of the run-and-tumble term for the case where θ=0 and ∇ϕ=∂ϕ/∂yŷ For this case,

(A9)

If I look at the change in θ over a small time period, τ, Δθ≈−(V₀τ/ϕ)(∂ϕ/∂y). Assuming that the bacterial velocity is v_b=V₀d and expanding to first order in Δθ, the continuity equation for the volume fraction 6 becomes

(A10)

In other words, the bacteria diffuse with diffusion coefficient , as is traditionally assumed 31. Therefore, this term will behave similarly to diffusion.

This model is similar to one developed by Alt 33 except that here it is assumed that the angular distribution is a constant normally distributed function, whereas Alt and Brosilow et al. 34 solve for the spatial and temporal evolution of the distribution. In addition, I assume here that during tumbling, the bacteria remain at the same spatial location.

Appendix B

Numerical implementation of the model equations

For all simulations, I set η_b=0. This allows me to use Eq. (16) to eliminate v_b,

(B1)

I also add Eqs. (6) to get a conservation equation that determines P:

(B2)

Equations (B1) are used in conjuncture with Eqs. (6) to solve for ϕ, the two components of v_f, θ, and P. However, I found that direct implementation of these equations into the general partial differential equation solver in COMSOL Multiphysics was unstable. I attributed this instability to the incompressibility condition (Eq. (B2)) and chose to use a method similar to artificial incompressibility to solve the system of equations 35. I therefore added time dependence to Eqs. (19) as

(B3)

where ɛ₁ and ɛ₂ are small. For this scheme, I need ɛ₂ to be less than ɛ₁ so that the fluid velocity always feels the steady-state, incompressibility condition (Eq. (B2)). Chorin 35 suggests that ɛ₂≤ɛ₁/10. Likewise, ɛ₁≪1 guarantees that the ϕ dynamics see the equilibrated fluid velocity. I found that for values of ɛ₁<0.01 and ɛ₂=ɛ₁/10 the solution was well converged.