Interesting results on the linear bioconvective instability of a suspension of gravitactic microorganisms have been calculated. The hydrodynamic stability is characterized by dimensionless parameters such as the bioconvection Rayleigh number R, the gyrotaxis number G, the motility of microorganisms d, and the wavenumber k of the perturbation. Analytical and numerical solutions are calculated. The analytical one is an asymptotic solution for small wavenumbers (and for any motility number) which agrees very well with the numerical solutions. Two numerical methods are used for the sake of comparison. They are a shooting method and a Galerkin method. Marginal curves of R against k for fixed values of d and G are presented along with curves corresponding to the variation of the critical values of Rc and kc. Moreover, those critical values are compared with the experimental data reported in the literature, where the gyrotactic algae Chlamydomonas nivalis is the suspended microorganism. It is shown that the agreement between the present theoretical results and the experiments is very good.
- hydrodynamic stability
- Galerkin method
Since many years ago, efforts in the experimental and theoretical investigation of the bioconvection phenomenon have been made. These efforts, which lead to the understanding of bioconvective instability, have produced novel and interesting applications. For example, Noever and Matsos  proposed a biosensor for monitoring the heavy metal Cadmium based in bioconvective patterns as redundant technique for analysis, a number of researchers [2, 3, 4, 5, 6] have been working on the control of bioconvection by applying electrical fields (as in galvanotaxis) to use it as a live micromechanical system to handle small objects immersing in suspensions, Itoh et al. [7, 8] use some ideas of bioconvection in a study for the motion control of microorganism groups like
The term bioconvection was first coined by Platt  as the spontaneous pattern formation in suspensions of swimming microorganisms. This phenomenon has some similarity with Rayleigh-Benard convection but originates solely from diffusion and the swimming of the organisms. Reviews about this topic have been published by Pedley and Kessler  and Hill and Pedley . Ideas and theories on cellular motility can be found in the book of Murase , and the effect of gravity on the behavior of microorganisms is widely explained in the book of Hader et al. . In 1975, Childress et al.  presented a model for bioconvection of purely gravitactic microorganisms and their results of a linear theoretical study, and later Harashima et al.  studied the nonlinear equations of this model. According to the model of Childress et al. , the critical wavenumber at the onset of the instability is always zero. In ordinary particles and colloidal suspensions, the internal degrees of freedom like the internal rotation or spin are important under some geometrical and physical conditions [21, 22]. The case of a suspension of microorganisms is not an exception. For this case, Pedley et al.  proposed a gyrotactic model for a suspension of infinite depth. Their model includes the displacement of the gravity from the geometric center in the organisms along their axis of symmetry. Hill et al.  performed an analysis of the linear instability of a suspension of gyrotactic microorganisms of finite depth using the model of Pedley et al. . Hill et al.  found finite wavenumbers at the onset of the instability, but agreement with the experiment was not good. Later, Pedley and Kessler  reported a model for suspensions of gyrotactic microorganisms where account was taken of randomness in the swimming direction of the cells. In a study of the linear instability of the system based on the model of Pedley and Kessler , Bees and Hill  found disagreement between their theoretical results and the experimental data reported by Bees and Hill . Several experimental investigations of bioconvection have been reported by Loeffer and Mefferd  and Fornshell , by Kessler  and Bees and Hill  who take into account the gyrotaxis, by Dombrowski et al.  and Tuval et al.  who take into account the oxitaxis, and more recently by Akiyama et al.  who observed a pattern alteration response characterized by a rapid decrease in the bioconvective patterns. Pattern formation has been observed in cultures of different microorganisms such as
More recently, investigations have been reported for a semi-dilute suspension of swimming microorganisms where cell–cell interactions are considered [34, 35, 36, 37, 38]. On the other hand, Kitsunezaki et al.  investigated the effect of oxygen and depth on bioconvective patterns in suspensions with high concentrations of
This chapter presents interesting results about the bioconvective linear stability of a suspension of swimming microorganisms. Use is made of the equations presented by Ghorai and Hill [48, 49] some years ago. In their approach, Ghorai and Hill [48, 49] used a different dimensionalization scale for the concentration microorganisms which gives distinct meaning to the basic state for the concentration of microorganisms and a bioconvective Rayleigh number defined in terms of the mean cell concentration. To the authors best knowledge, those equations along with the change in the basic state and Rayleigh number definitions have not been used to determine the linear bioconvective instability in an infinite horizontal fluid layer and to compare the results with experiment. These results were obtained by means of both numerical and analytical techniques. The critical values of the Rayleigh number and the wavenumber , for fixed values of the gyrotaxis number
The chapter is organized as follows. The governing equations and boundary conditions [48, 49] as well as the basic state can be found in Section 2. Nondimensionalization and linearization of the system of equations is outlined in Section 3. In Section 4, use is made of an asymptotic expansion [50, 51, 52, 53] method and a Galerkin method  to find limiting cases and predict critical values of
2. Equations of motion
We consider an infinite horizontal layer of a suspension of gyrotactic microorganisms. The fluid layer is bounded at . The fluid where the cellular microorganisms swim is water with density . Each cell has a volume and density , where . The suspension is considered dilute and incompressible. Density fluctuations in the suspension are small enough such that the Boussinesq approximation is valid and the corresponding governing equations are
where is the time,
where is the cell swimming speed,
where represents the average concentration of organisms. Eq. (7) is the same basic state as presented by Ghorai and Hill [48, 49] whose linear stability will be investigated. It differs from that of Childress et al.  and Hill et al.  by the coefficient
3. Linear stability
We make the governing Eqs. (1–3) nondimensional by scaling all lengths with
where the nondimensional quantity is the ratio of swimming speed of microorganisms and their representative mass diffusion velocity. Here,
are the Schmidt and bioconvection Rayleigh numbers, respectively. Pedley and Kessler  give a definition of the vector
where the subscript denotes the horizontal component, is the cell eccentricity, and is the nondimensional form of the gyrotactic orientation parameter . Finally after substitution of
with boundary conditions
where the superscripts have been deleted. Notice that the adimensionalization of the equations is different from that of Hill et al. . Here, an application of a more general asymptotic analysis for any magnitude of
By elimination of the pressure from Eqs. (13–15), it is possible to obtain a coupled system of two equations, for
where is the wavenumber of the disturbance and is the growth rate. The wavenumber is scaled as corresponding to a nondimensional wavelength . Thus, the governing equations become
subject to the boundary conditions
where . The variables of the above problem can be changed in order to simplify the analysis. The change of dependent variable is
subject to the new boundary conditions
In this form, the equations are very similar to those of the well-known problem of thermal convection in an infinite horizontal fluid layer between nonconducting boundaries [50, 51, 52, 53, 56]. The familiar fixed heat flux boundary condition is the main characteristic of those thermal convection problems and is analogous to that presented in Eq. (22). The equations derived by Childress et al.  can also be analyzed from the present view point of this change of variable. In the theory of thermal convection as in that of Childress et al. , a zero critical wavenumber is found as a result of the fixed flux boundary condition. In more recent models, which include the effects of gyrotaxis, the similarity with the thermal convection problem is not valid unless .
4. Asymptotic analysis
In this section, the eigenvalue problem stated in the system of Eqs. (13–15) with boundary condition Eq. (16) is investigated by means of two analytic methods. The magnitude of the marginal value of
4.1 Asymptotic analysis
We conducted a general asymptotic analysis in comparison with those used before [19, 24, 26] which included the restrictions of the limits for shallow layers and for deep layers along with different restrictions for
where . We also consider no restrictions for and and rescale the wavenumber as . Thus, after substitution of expansions Eqs. (23–26) and the mentioned scalings in Eqs. (20) and (21) with boundary condition Eq. (22), we obtain the following systems of equations at different orders.
The systems of equations at order and higher are inhomogeneous and must satisfy their corresponding solvability conditions allowing to compute the Rayleigh number
Thus, the solvability conditions at and are, respectively,
The solutions of the system of equations at leading order are
where the function can be obtained from the authors upon request. For convenience, the solution has been normalized to 1. The next step is to evaluate the solvability condition Eq. (36) at and obtain an expression for
The constant is large and can be obtained from the authors upon request. At order similar steps as those for solving the system of equations at are followed to find and . Then, algebraically is
After substitution of , , , and into Eq. (30), the velocity can be calculated subject to its corresponding boundary condition Eq. (32). Because of the term appearing in the system of equations at , the expression of is very large and complicated and will not be given here. The evaluation of the solvability condition at order given in Eq. (37) yields
The growth rate can now be obtained by substitution of and into the expansion for given in Eq. (26). However is omitted to save space but can be obtained from the authors upon request. Finally, use is made of the expansion Eq. (25) for .
Now, the transition from stability to instability via a stationary state is investigated by setting . Thus, the corresponding value of for the marginal state is
where some simplifications have been made with the use of obtained from Eq. (39). The functions and and the coefficients to appearing in the above expressions are functions of and can be obtained from the authors upon request. The result for is
From the expression for the Rayleigh number given in Eq. (42), it is possible to calculate the limit for . In this case, we consider and to be of the same order in such a way that and . Then, under these assumptions, the approximation of is
Here we point out that in the present chapter, our definition of the Rayleigh number differs from that defined by Hill et al. . If our approximation given in Eq. (44) is multiplied by , the same approximation by Hill et al.  is obtained. Moreover, if this becomes that given by Childress et al. . In the more general expression of for a small wavenumber approximation, Eq. (42) has a special characteristic due to its dependence on the square of the wavenumber . The first coefficient at zeroth order in corresponds to , and that at the second order in is . Even though in the experiments on bioconvection , only finite critical wavenumbers have found these coefficients are very useful. For example, it can be shown from that if and
then . This corresponds to a stable stratification, which is not the case here. The second coefficient in Eq. (42) is a very important result, because it provides information about the shape of the marginal curve with respect the critical wavenumber. That information can be obtained by making zero the coefficient and calculating the following critical value of the gyrotaxis number
From this equation, two admissible cases are possible when Eq. (45) is satisfied. First, for fixed values of , , and , the marginal curve starts at and then decreases monotonically. However, according to the numerical analysis presented below, the marginal curves in fact first decrease and then start to grow monotonically after a minimum is attained, at the critical wavenumber. In the second case, for fixed values of , , and , the marginal curves start at and then grow monotonically. Here, the critical wavenumber is always zero. The importance of these results is that the magnitude of agrees very well with the results of the marginal curves found in the numerical analysis given below. This critical value determines the magnitude for which the curves have a finite critical wave number () or a zero critical wavenumber (). It is of interest to note that some of the magnitudes of the gyrotaxis parameter calculated from the data in the literature are very near but above of . This is the reason why some of the curves found from numerical analysis are almost flat in a range of wavenumbers near to zero. Because the experimental critical wavenumbers found for gyrotactic bioconvection are always finite, we conclude that is important to point out where to find the theoretical limitations of the model.
4.2 Analytic Galerkin method
Here use is made of the analytical Galerkin method to study the eigenvalue problem of Eqs. (17)–(18) with the boundary condition Eq. (19). This method has been used before by Pellew and Soutwell , Chandrasekhar , and Gershuni and Zhukovitskii . Even though this is an approximate method, it has a very high precision. The advantage of the method is that it is possible to obtain an explicit expression of the Rayleigh number
Briefly, the method consists in assuming a trial function which satisfies the boundary conditions for each of the dependent variables. Let that variable be which after substitution in one of the equations of the problem allows for the exact solution of the other variables, let us say
In this way, the proposed expansions of and are
which is subjected to the conditions
The solution is
where to are constants of integration which can be found by evaluating in the boundary condition Eq. (49).
Next, Eq. (18) is multiplied by and is integrated in the range to , to get
This determinant, calculated with the help of the software Maple, is the solvability condition from which the eigenvalue is calculated. The resulting algebraic expression of the integrals in this equation is very complex and will not be presented here. However, the first approximation of , corresponding to the element (0,0) of the matrix in the determinant Eq. (52), is
where and are functions of the wavenumber
5. Numerical computations by a shooting method
Here, the shooting method  is used to solve the eigenvalue problem posed by the system of Eqs. (20) and (21) subjected to the boundary condition Eq. (22). Curves of marginal stability in the plane were calculated for fixed values of the parameters , , and . Notice that very good agreement was always found among the values of the of the asymptotic analysis (in the limit ), of the Galerkin method, and of the numerical computations. Calculations were made in two ways. First, the parameters , , and were varied in order to obtain a representative set of marginal curves for the problem of bioconvection. Second, experimental data were also used to fix the values of , , and and used to find theoretical values of and that could be compared with their corresponding experimental values. Here, in particular, a selection is made of , which corresponds to the flagellated alga
In the curves shown in Figure 1a–b, the critical values of the gyrotaxis parameter are , respectively. As mentioned above for , the critical wavenumber is finite, and for the critical wavenumber is always zero. The combined effects of the velocity of the swimming of microorganisms, , and that of gyrotaxis, , change the location of the critical wavenumber. Note also that for fixed , when increases, the system becomes more unstable. From the Figure 1a and b, it can be seen that the most unstable case corresponds to that for and where and . This may be understood by the fact that the accumulation of microorganisms near to the top of a shallow layer is faster than in a deeper one. This is due to the important role that the mass diffusion of microorganism and the depth of suspension play on the instability of the system. The value of in the limit of can also be calculated from Eq. (46) by means of an asymptotic analysis. That is,
Here, some theoretical curves are presented of which some have a very good agreement and others a reasonable agreement with the experiments 2, 4, 9, 10, 13, 16, 20, 24, 26, 27, 28, 29, 31, and 35, performed by Bees and Hill .
The values for the motility and the gyrotaxis parameter used in Figures 2 and 3 were calculated based on experimental data by Bees and Hill [26, 27], which in here are presented in Tables 2 and 3 of the following section. In order to observe in detail the position of the critical point in Figures 2 and 3, a local magnification is included.
Here, a comparison is done of our theoretical results of with the theoretical ones presented by Bees and Hill  in their Table VI. According to Bees and Hill , experiments 2 and 23 in their Table V have of comparable order with those in their Table VI. In our Table 1, we reproduce the comparison made by Bees and Hill  of their own theoretical and experimental results of their Table V, and we added the corresponding error in percent of the wavenumbers and Rayleigh numbers, respectively. Note that the value corresponds to flagellated microorganisms such as
|Experimental results||Theoretical predictions||Error (%)|
For the sake of comparison of our theoretical results with those of the experiments, Table 1 shows the percent of error calculated by taking the difference of the experimental and theoretical values and then dividing by the smallest one. In Table 1, the more realistic value for
6. Comparison with experiments
In this section a comparison is done of our theoretical results of and with the corresponding experimental values obtained by Bees and Hill . Here use is made of the results of the 39 experiments shown in Table I of Bees and Hill . Besides, the more realistic value of the parameter , corresponding to the flagellated algae
In Table 3, the values of , , and resulting from the experimental data are presented. Note in Table 2 that the experimental results of the cell swimming speed and of the cell diffusivity are given inside a range of values. In this case, a particular value inside the range has to be selected. The swimming speed used here is 63. The decision is based on the results obtained by Hill and Hader , Pedley and Kessler , and Bees and Hill . The value of the cell diffusivity was decided to be that corresponding to an average over the range given in Table 2, that is, .
|Cell volume||5 × 10−10 cm3|
|Acceleration due to gravity||103 cms−2|
|Cell diffusivity||5 × 10−5–5 × 10−4 cm2s−1|
|Fluid density||1 gcm−3|
|Cell density||1.05 gcm−3|
|Kinematic viscosity||10−2 cms−2|
|Cell swimming speed||0–2 × 10−2 cms−1|
|Dimensional gyrotaxis parameter||3.4 s|
|Including flagella||6.3 s|
Very recent experimental measurements on the diffusivity for different microorganisms like the biflagellated alga
where the constants and (see  for more details). is a direction correlation time which equals in the nonflagellated case and in the flagellated case. The data corresponding to the suspension depth and the average cell concentration of microorganisms of each experiment (see  for more details) have not been reported in Table 3. Only the parameters , , , and are presented in that table. It is also found that the value of the of each experiment is greater (but sometimes near) than their corresponding critical value of Eq. (46). Under these conditions, all the critical wavenumbers have to be .
|NE||Error ||Error |
By using the data of our Table 2 and Table I of Bees and Hill , the experimental values for
Some numerical results agree very well with experiments, as can be seen in the experiments 4, 10, 12, 13, 20, and 35 of Table 3. Others are good, such as the results of experiments 9, 16, 24, 26, 27, 28, 29, and 31. With respect to the other data in Table 3, it might be possible that if the experimental measurements are improved, the agreement with theory will be better. The results given here show that the approximate and numerical solutions of the system of governing equations presented in this paper may bring a light to the solution of many other problems of bioconvection.
The governing equations of bioconvection were used to investigate the problem of an infinite horizontal microorganism suspension fluid layer. The theoretical predictions of the critical wavenumber and Rayleigh number were compared with their experimental counterparts . Very good, good, and fair agreements were found. But in general, we may say that our numerical results improve by far those obtained by Bees and Hill .
With the asymptotic analysis for , it was possible to calculate a Rayleigh number not reported before without any restrictions on the magnitudes of
However, it is clear from the experimental results that the critical wavenumbers are finite and large and that the former case is not physical. Therefore, this also defines the limit of validity of the theory. Note that it agrees very well with numerical analysis.
An analytic Galerkin method was also used to obtain a general expression of
Numerical results have shown that the system becomes more unstable when the layers are shallow. The physical interpretation of such situation is that the accumulation of microorganisms near the top of the layer in the shallow case is faster than in the deeper case, due to the smaller depth of suspension
Finally, we would like to point out that it is our hope that the results presented in this chapter may stimulate researchers to make more new and precise experiments on bioconvection in the near future.
dimensional gyrotactic parameter, s wavenumber cell diffusivity, cm2s−1 motility of microorganisms dimensionless gyrotactic parameter acceleration due to gravity, cms−2 layer depth, cm flux density of organisms average cell concentration concentration of microorganisms pressure Rayleigh number Schmidt number time cell swimming speed, cms−1 fluid velocity Cartesian coordinates cell eccentricity viscosity, gcm−1s−1 kinematic viscosity, cm2s−1 water density, gcm−3 cell volume, cm−3 result obtained by Bees and Hill  critical value experimental result theoretical result
dimensional gyrotactic parameter, s
cell diffusivity, cm2s−1
motility of microorganisms
dimensionless gyrotactic parameter
acceleration due to gravity, cms−2
layer depth, cm
flux density of organisms
average cell concentration
concentration of microorganisms
cell swimming speed, cms−1
kinematic viscosity, cm2s−1
water density, gcm−3
cell volume, cm−3
result obtained by Bees and Hill 
The authors would like to thank Alberto López, Alejandro Pompa, Cain González, Raúl Reyes, Ma. Teresa Vázquez, and Oralia Jiménez for technical support. I. Pérez Reyes would like to thank the Programa de Mejoramiento del Profesorado (PROMEP).