Environmental Transport of Gyrotactic Microorganisms in an Open‐Channel Flow

The rich and complex phenomena of the transport of active particles like microorganisms in shear flows are of great significance to various biological and environmental applications. Recent studies have shown that the motility and gyrotaxis of algae could greatly influence their transport in waters. However, little attention has been paid to the initial and transient transport regime when the classical Taylor dispersion model is not applicable. To tackle this problem, we resort to Gill's generalized dispersion model for passive particles like solute, which has the potential for accurately describing the entire transport process. For the first time, we extend Gill's model to the active particles, and the effects of swimming, gyrotaxis, and flow shear on the microorganism dispersion in an open‐channel flow have been thoroughly investigated. We first theoretically solve the transient drift and dispersion coefficients, based on which we obtain analytical solutions for concentration distributions of microorganisms, and further validate these results by numerically solving the governing equation. We find that when there is no flow, the longitudinal dispersion of microorganisms can be weakened by the gravitactic accumulation in the vicinity of water surface, while enhanced by a stronger swimming ability of the microorganisms. The effect of the flow shear does not affect the form of the asymptotic concentration distribution, but can greatly enhance the transient drift velocity and the dispersivity. We further analyze the effect of turbulence on microorganisms' dispersion by combing the direct numerical simulation and the random walk simulation. The increase of turbulence is shown to decrease the vertical non‐uniformity of the concentration distribution, as well as the relative contribution of active behavior to both the drift and Taylor dispersivity during transport.

However, both the PK model and the GTD method are subject to restrictions resulting from the assumption of separable orientational and translational dynamics. For example, the diffusion tensor obtained by the PK model is only an acceptable approximation in weak shear (Pedley & Kessler, 1990), while that obtained by the GTD method is only valid for flows with shear weakly varies spatially (Bearon et al., 2011(Bearon et al., , 2012. Since these conditions are normally not satisfied in typical confined shear flows, the effective dispersivity obtained for the dispersion of gyrotactic microorganisms, based on the approximate diffusion tensor, could deviate from the exact result. Recently, W. Jiang and Chen (2019) proposed a rigorous extension of the original GTD theory (Frankel & Brenner, 1989) for the effective one-dimensional longitudinal dispersion of micro-swimmers: they assigned all bounded coordinate(s) as the local space and all unbounded coordinate(s) as the global space, which was previously suggested in the paradigmatic format of the GTD theory (Brenner & Edwards, 1993). The GTD theory was then applied to the asymptotic dispersion in 2-D channel flows (W. Jiang & Chen, 2019) and pipe flows (W. Jiang & Chen, 2020;Wang et al., 2022a), and the results were validated against the random walk simulations . For the dispersion of active particles in bounded flow, the previous research mainly focused on the long-time asymptotic characteristics. Recently, W. Jiang and Chen (2021) developed a semi-analytical method to calculate the temporal concentration moments, by introducing the biorthogonal expansion method.
There is a growing interest on investigating the transport of gyrotactic microorganisms in flows with free-surface (Guan et al., 2021;Lovecchio et al., 2017;Mashayekhpour et al., 2019;Yang et al., 2020Yang et al., , 2021Zeng et al., 2019). Nonetheless, these works either used the random walk numerical simulation or the further simplified version of the PK model. It is thus interesting to see if theoretical analysis can be performed to more accurately describe the transport of gyrotactic microorganisms in shear flow with a free-water surface.
To address the transient evolution of dispersion instead of just focusing on the asymptotic process, one can resort to the classic generalized dispersion model proposed by Gill (1967), initially devised and currently applied for the transport of passive particles. The second-order generalized dispersion model is in the same form as the Taylor dispersion model, but with time-dependent coefficients for the advective and diffusive terms, respectively, which can capture the exact information regarding the transient mean velocity and the expansion of the concentration cloud. Although higher-order terms can be incorporated into the model (Gill & Sankarasubramanian, 1970) to capture more statistical information, the second-order generalized dispersion model is most widely applied in 10.1029/2022WR033229 3 of 26 different configurations (Guo et al., 2022;W. Q. Jiang & Chen, 2018;Rana & Murthy, 2017;Sadeghi et al., 2020;Wu et al., 2015) due to the feasibility in analytically deducing the required time-dependent coefficients.
The rest of this paper is structured as follows. For the transport problem of gyrotactic microorganisms formulated in Section 2, we extend Gill's generalized dispersion model to the transport of active particles in Section 3. The corresponding governing equations are solved by series expansion on the basis of the longitudinal derivatives of cross-sectional-mean concentration. We validate the obtained analytical solutions by numerically solving the governing equation using random walk simulations in Section 4. Aiming at proposing a general framework instead of analyzing for a specific application, we investigate the effects of swimming, gyrotaxis, and shear flow on the transient dispersion process in an laminar open-channel flow, as shown in Section 5. The temporal dispersion characteristics such as drift and the effective dispersion coefficient are solved, based on which the evolution of both the vertical-mean and the two-dimensional concentration distribution are obtained and discussed. A primitive discussion on the effect of turbulent flow is provided in Section 6. Finally, Section 7 gives some conclusions and prospects.

Governing Equations
We consider an extremely basic transport problem to analytically study the transient dispersion process. As shown in Figure 1, a dilute suspension of gyrotactic microorganisms is dispersing in a two-dimensional laminar open-channel flow. The transport process is affected not only by the translational diffusion and the advection of the shear flow, but also by the rotational diffusion and swimming of the microorganisms. Therefore, apart from the position coordinates (vertical coordinate z* and longitudinal coordinate x*, where the superscript * indicates dimensional variables), an orientation coordinate θ indicating the swimming direction should also be taken into consideration.
For simplicity, only the laminar background flow is taken into analytical analysis. Although the turbulent flow is also of interest to environmental applications, it would greatly increase the complexity. In Section 6, we have discussed the turbulent effect by numerical simulations and have found that the laminar case can capture the key evolution of the active transport characteristics and thus is appropriate for the following theoretical study. Additionally, because of the dilute suspension, the cell-cell and cell-fluid interactions can be ignored in the transport model. It is also assumed that the swimming velocity of microorganisms is a constant, denoted by * .
Based on the above assumptions and conditions, the governing equation for the transport can be adopted as (Wang et al., 2022b) where t is time and P(x, z, θ, t) is the probability density function (pdf) of an individual of the microorganism found in the position and orientation space {x, z, θ, t}. In Equation 1 we have already used the following dimensionless variables and parameters 10.1029/2022WR033229 4 of 26 the rotational diffusivity. The translational diffusion is negligible compared to the swimming-induced diffusion (Bearon et al., 2011), and thus D t = 0.
Note that the definition of the dimensionless longitudinal coordinate x suggests that the coordinate system is transformed into a reference frame moving with the mean flow speed * ≜ 1 where H* is the flow height. The laminar open-channel flow profile is a parabolic function of z, with a no-slip boundary at the bottom and a free-surface condition. The dimensionless flow velocity minus its mean is Note that in Equation 1, ̇ is the change rate of the swimming direction. Apart from the change induced by the ambient shear flow, the swimming direction is biased by the gravitational torque, which is the typical characteristic of gyrotactic microorganisms. ̇ can be calculated by the Jeffery's equation (Hinch & Leal, 1972;Jeffery, 1922) = 1 2 (−1 + 0 cos 2 ) + cos (5) here α 0 is the shape factor of the microorganism ranging from 0 (spherical) to 1 (rod-like, infinitely long); λ is the dimensionless bias parameter with B* representing the gyrotactic reorientation time (Pedley & Kessler, 1987).

Initial and Boundary Conditions
Next, we specify the initial condition and boundary conditions for the governing Equation 1. For simplicity, the microorganisms are released at the middle of the channel, and the initial distribution of the swimming direction is uniform. Thus, the initial condition of P is For the boundary conditions, it is assumed that the collision between the microorganisms and the bottom of the channel is completely elastic (Bearon et al., 2011). Because microorganisms cannot swim out of the water, we also make an idealized elastic collision hypothesis at the free surface in order to facilitate calculation (Mashayekhpour et al., 2019). Then, the boundary conditions for P (x, z, θ, t) can be written as which ensure the zero-total-flux condition at the boundaries: Thus, the total amount of microorganisms is conserved during the transport. This kind of boundary condition is called a reflective condition (Bearon et al., 2011;Ezhilan & Saintillan, 2015).
We remark that the real particle-wall collision is much more complicated than the above elastic assumption. In fact, the reflective boundary condition cannot account for some distinctive features of active suspensions in confined environments, for example, wall accumulation. Ezhilan and Saintillan (2015) have demonstrated the boundary-layer-like accumulation effect using a Robin boundary (Elgeti & Gompper, 2013;Enculescu & Stark, 2011;W. Jiang & Chen, 2019) with a non-zero translation diffusion coefficient. H. Chen and Thiffeault (2021) have further considered the geometrical constrain by microorganism shape onto the boundary condition. Moreover, there are many other models proposed to address the accumulation by hydrodynamic interactions (Berke et al., 2008;Spagnolie & Lauga, 2012), steric inter-molecular forces like the van der Waals force (Contino et al., 2015;Costanzo et al., 2012;Li et al., 2008), and the change of swimming behaviors . Imposing these complicated models in analytical studies is challenging. Fortunately, in the current setup, the surface accumulation by gravitaxis is stronger than the wall accumulation, as the vertical concentration distribution is nearly exponential. Besides, the translational diffusion is often ignored in many studies, so we focus on the more typical reflective condition, which is also widely accepted in numerical studies.
For the angle of the swimming direction, the periodic conditions are automatically satisfied

Gill's Series Expansion for Swimming Microorganisms
The generalized dispersion model by Gill (1967) provides one of the most classic analytical methods to study the dispersion of passive particles. Here we extend the model for the transport of swimming gyrotactic microorganisms. The essential idea of the generalized dispersion model is to separate variables and expand the concentration distribution C into a series of the longitudinal derivatives of the cross-sectional mean concentration ̄ (Gill & Sankarasubramanian, 1970, Equation 5) where f n are the expansion coefficient functions.
For swimming microorganisms, we need to include the variation of the swimming directions. An analogous series expansion for P (x, z, θ, t) can be introduced as where the overall mean of P is defined by the angle brackets as which describes the evolution of the longitudinal distribution for the vertical-mean concentration of the microorganisms. f n are functions of z, θ, and t, and by taking the overall mean as defined by Equation 12, we find that Substituting the pdf expansion Equation 12 into the governing Equation 1 and performing the overall mean operation, we obtain the important Taylor-Gill expansion equation: 10.1029/2022WR033229 6 of 26 where the key transport coefficients can be expressed as and f −1 = f −2 = K 0 = 0.
Assigning n = 1 for Equation 16, K 1 is found as which represents the drift speed of the overall microorganism cloud.
Similarly, assigning n = 2 for Equation 16, K 2 is found as which represents the Taylor dispersivity (or the effective diffusion coefficient).
Equation 15 is the general dispersion model for swimming microorganisms, in the same form as that for passive particles (Gill & Sankarasubramanian, 1970, Equation 8). As can be seen in Equations 17 and 18, both K 1 and K 2 are time-dependent compared with the constant coefficients in Taylor dispersion model, demonstrating the ability to describe the transient transport process. Before determining Equation 15, we need first to solve the transport coefficients f n and K n .

Solutions of Transport Coefficients
The governing equation for f n (z, θ, t) can be obtained by substituting the series expansions Equations 12 and 15 into Equation 1. Collecting the expanded results in terms of the longitudinal derivatives, we have The initial condition for f n can be given according to Equation 7, as The boundary conditions for f n are in the same form as Equation 8, (1, , ) = − (1, − , ).
(21d) f n and K n can then be solved by the Galerkin method, which is similar to that used in previous papers (W. Wang et al., 2022b). Below is a summary of the solution procedure, and more details can be found in Section 3.2 of the work by W. Jiang and Chen (2021).
The basic idea is to express f n with a series expansion as where { } ∞ =1 is an orthogonal function basis and a ni is the expansion coefficient. According to the boundary conditions Equation 21 for f n , the reflective basis (W. Jiang & Chen, 2021, Equation 3.26) can be used, which is comprised of where i = 1, 2, … and j = 1, 2, …. We adopt the first 40 functions of z and the first 10 functions of θ as a truncation of the basis functions for the summation of ψ in Equation 22. Substitute Equation 22 into the governing Equation 19 and use the orthogonal relationship, we can obtain the ordinary differential equations of a ni , which is a system with constant coefficients and thus can be solved easily by matrix diagonalization with corresponding eigenvalues. After solving f n , one can obtain K n by Equation 16 immediately.

Solutions of Generalized Dispersion Model
Here we solve the Taylor-Gill expansion Equation 15 at order two, which is also the Taylor dispersion model for an approximation of mean concentration 〈P〉. With K n (t) and f n (r, t) solved in the above section, we substitute them into Equation 24, and the corresponding solution for such a classic advection-diffusion equation can be given based on the fundamental solution of the diffusion equation: where ω n is the integral of K n : which is a Gaussian distribution.
To obtain an approximation for the pdf containing the information on vertical positions and swimming directions of microorganisms, we resort to the generalized dispersion model (11) at order 1: The two-dimensional concentration distribution is obtained by considering swimming microorganisms toward all directions: where 〈⋅〉 θ denotes the integration over all swimming directions. Furthermore, one can obtain the vertical concentration distribution of microorganisms by integrating out the dependencies on both x and θ coordinates:

Comparison With Numerical Simulations
To verify the obtained analytical results based on Gill's dispersion model, we perform the widely used random walk simulation to numerically solve the governing Equation 1. The corresponding stochastic differential equations for Equation 1 are where W x , W z , and W θ are independent standard Brownian motions.
In the random walk simulation, we adopt a small enough time step of Δt = 2 × 10 −4 , and a total number of 10 5 microorganisms are considered to ensure a relatively good convergence of the simulated results. As a demonstration of the performance of the obtained analytical solutions, here we compare both analytical and numerical results for the case with the specific set of parameters: Pe f = 0.2, Pe s = 0.5, and λ = 1, under which circumstance the effects of the flow shear, the microorganism swimming, and the gyrotaxis are all present.
The most important message we want to convey by Figure 2 is that the transport coefficient K n obtained in our approach contains the exact information regarding the statistics for the spatial-temporal evolution of the microorganism cloud. Specifically, for Gill's dispersion model by the second order as used in this study, −K 1 represents the downstream velocity of the centroid of the cloud, and K 2 indicates the ability for the streamwise scattering of the microorganisms, both of which show a perfect match between the analytical and numerical results as illustrated in Figure 2. We note that the classic Taylor dispersion model is a one-dimensional advection-diffusion equation with the coefficients K 1 and K 2 being constants, equivalent to the currently used model for relatively large times (roughly t > 6 for Figure 2), which is not applicable to the earlier transient stage for the transport of microorganisms.
Rather than the exact one-dimensional (streamwise) statistical information contained in K 1 and K 2 , we further explore the spatial concentration distribution itself for the microorganism cloud. Figure 3 presents the two-dimensional distributions for both analytical and numerical results at different times. It is seen that at the early stage the analytical results largely capture the "real distribution" as represented by the numerical simulation, but differences between the two are also obvious. As time goes on the two distributions are more and more alike, which become almost identical for Figures 3c and 3f. These observed differences are expected since we considered only the first two orders of the expansion for the mean concentration in Equation 15, and only the first order of the expansion for the transverse concentration distribution in Equation 27.
For a more detailed understanding for the concentration distribution, we continue to display the streamwise distribution of the vertical-mean concentration in Figure 4. Here we can see that typical non-Gaussian characteristics exist for distributions at the early stage of the transport. Although the analytically predicted Gaussian curve fails to describe the nearly rectangular distribution by the numerical simulation in Figure 4a, it correctly captures the streamwise scale of the concentration cloud. This is due to the correctly predicted K 2 at the early time, based on which the streamwise scale of the cloud can be approximated by the standard deviation as √ 2 2 . These non-Gaussian characteristics of course are related to higher-order terms in the expansion for the mean concentration (Equation 15).
Compared to the relatively poor performance for the mean concentration distribution prediction, in Figure 5, we find in surprise that the analytical solution for the vertical concentration distribution agrees perfectly with the numerical simulation, even for the early stage of the transport (Figure 5a). This is important for later in this paper analyzing the gyrotactic trapping of the microorganisms at a certain depth of the flow (Figure 5a as an example) using the analytical solution.
Overall, by comparisons with numerical simulations, we demonstrate that the obtained analytical solutions contain accurate statistical information regarding the first two orders of concentration moments. They describe with reasonable accuracy for the two-dimensional concentration distribution at the initial transient stage of the microorganism transport, mainly due to the truncation error for the vertical-mean concentration. Instead, the vertical concentration distribution can be well predicted by the analytical solution.

Effects of Gyrotaxis, Swimming, and Flow Shear
After demonstrating the capability of the obtained analytical solutions in describing the microorganism concentration distributions in laminar open channel flows, in this section, we discuss separately the key factors affecting their transport process, including the effects of the gyrotaxis, the swimming, and the flow shear, respectively, based on the Taylor-Gill dispersion model extended in this work.

Influence of Gyrotaxis
The bias parameter λ reflects the strength of gyrotaxis, which describes the upward swimming tendency of the microorganisms, and further determines the degree of gravitactic accumulation near the water surface in the open-channel flow. In order to analyze the effect of gyrotaxis on the dispersion process, we consider the microorganisms as spherical  particles with different bias parameters, given a constant swimming speed and zero velocity for the ambient fluid flow. Specifically, we take four different values of 0, 0.5, 1, and 2 for λ, and keep Pe s = 0.5 and Pe f = 0 unchanged.

Dispersivity
As clearly indicated in Equation 17, the drift velocity K 1 includes two parts: the advective part contributed by the fluid flow, and the swimming part by the active motions of the microorganism itself. With no motions of the ambient fluid flow, f 0 is symmetric about θ = π/2; thus we can see immediately that the overall drift speed of the microorganism cloud is 0 regardless of the value of Pe s , which indicates the swimming ability of the microorganisms.
The effective diffusion coefficient K 2 characterizes the longitudinal scattering of the cloud of microorganisms. As shown in Figure 6, K 2 eventually reaches a steady value at t ∼ 5, demonstrating that the asymptotic dispersion regime is attained. For the influence of gyrotaxis, we see that as λ gets larger, K 2 decreases for both the transient and steady values. It is understood that if the gyrotactic bias is larger, there are more microorganisms swimming toward the water surface while less toward other directions, thus limiting the overall longitudinal dispersion of microorganisms as reflected by a smaller K 2 .

Concentration Distribution
As shown in Figure 7, other than the special case of λ = 0 which is more like the diffusion of passive particles leading to a uniform vertical concentration distribution, microorganisms gradually migrate to the water surface due to gyrotaxis after being released at the middle of the channel. And as the gyrotactic bias λ becomes larger, more microorganisms get concentrated close to the vicinity of the water surface at the same time. Now we focus on the longitudinal distribution of vertical-mean concentration. As shown in Figure 8, the mean concentration is symmetrically distributed with respect to x = 0. As time increases, microorganisms diffuse to both sides of x = 0: the maximum of the mean concentration around x = 0 decreases, and the longitudinal distribution becomes flat. As reflected by a smaller K 2 under larger λ, the diffusion of the vertical-mean concentration is slower with a larger λ, resulting in a higher maximum and smaller longitudinal extent of the cloud at the same time. Again, this is because if the gyrotactic bias is larger, more microorganisms tend to swim toward the water surface, and therefore more microorganisms are found in the vicinity of the release point, corresponding to a larger concentration near x = 0.
We further discuss the evolution of the vertical concentration distribution C v , which is the integral of microorganisms concentration with respect to both x and θ, describing the distribution of microorganisms at each streamline. A really interesting phenomenon as observed in Figure 9b is that although microorganisms are initially released at z = 0.5, the concentration around that position becomes almost the lowest at t = 2 (e.g., a U-shape of the distribution in Figure 9b for λ = 0). This is very different from the diffusion of passive particles, for which case we always observe the maximum around the release position, which is also the symmetrical center of the Gaussian distribution. This can be explained by the swimming of the microorganisms: once they are released at z = 0.5, they immediately swim away (equally probable  up or down for the case of λ = 0), leading to a minimum of the concentration at the vicinity of the release position; as they approach the upper and lower boundaries they are subject to the reflective boundary conditions 8, which overall force them swimming back and smear out the concentration non-uniformity.
On the other hand, the asymmetry of the vertical distribution increases as the increase of λ (Figure 9a), due to the overall upward motions of the microorganisms, further resulting in the concentration higher close to the surface while lower around the bottom, which is a direct reflection of the effect of gravitactic accumulation. For the larger time of t = 10, the asymptotic vertical concentration distributions are displayed in Figure 9d, illustrating the vanishing of the temporally trapped microorganisms around the release position. It is obvious that the non-gyrotactic microorganisms (λ = 0) will eventually form a vertically uniform concentration distribution. Note that this uniform distribution is due to the reflective boundary conditions 8, which cannot account for the wall accumulation effect, as discussed in Section 2.2 and previous studies W. Jiang

Influence of Swimming
To analyze the influence of swimming on the dispersion process, spherical microorganisms with different swimming abilities (quantified by Pe s ) are released instantaneously in the middle of the cross-section of the channel. Specifically, the swimming Péclet numbers Pe s varies and the gyrotactic bias parameters are fixed (λ = 1). Pe s = 0 is equivalent to the case of passive particles. To understand the influence of swimming, we continue to impose no background shear flow (Pe f = 0).

Dispersivity
Similar to Section 5.1.1, with no background flow, f 0 is symmetric about θ = π/2, thus the overall drift speed of microorganisms is always 0.
Figure 10 demonstrates another important difference in the transport characteristics between active and passive particles. For the limiting case, if the microorganism cannot swim at all (Pe s = 0), there will be no difference in the pattern of scattering for the active and passive particles (Figure 10b), because their movements are only subjected to the (translational) Brownian motion and can be described by the pure diffusion process (constant K 2 all the time). The only difference in this case is the magnitude of the constant (K 2 ): for salt in still water, K 2 is the molecular diffusivity which is on the order of d − 8; for the microorganism it can be much smaller due to the relatively large cell compared to the salt molecule. For example, we can estimate that D t ≈ 6 × 10 −13 by the Stokes-Einstein relation for a spherical particle with diameter 10 μm in a 1-m-high channel at 20°C.
D t for the microorganism is so small that it can be considered negligible even compared to the diffusion of salt (if fact, we set D t = 0 in the numerical simulation and analytical solution). However, once the microorganism starts to swim, everything is different. First, the swimming is not a random motion, thus the scattering of the microorganism cloud is not a diffusion process until a certain timescale when K 2 becomes a constant (approximately by the order of 10s according to Figure 10a). Strictly speaking, the effective diffusion of active microorganisms originates from the combination of translational and rotational diffusion, and is affected by the presence of two boundaries. Second, the magnitude of K 2 increases greatly due to the swimming of microorganisms, here by approximately nine orders of magnitude for Pe s = 0.1 compared to that of Pe s = 0. And K 2 can still increase as the increase of Pe s , which scales approximately 2 ∼ 2 .

Concentration Distribution
As shown in Figure 11, as time increases, microorganisms migrate toward the water surface. After the same time, the aggregation of microorganisms with larger Pe s near the water surface is less obvious.
For longitudinal distribution, the vertical-mean concentration is axisymmetric regarding x = 0 where the microorganisms are released. As shown in Figure 12, with the increase of time, microorganisms swim and diffuse to both sides of x = 0; the longitudinal distribution of vertical-mean concentration becomes flat and the vertical-mean concentration near x = 0 decreases. Increasing Pe s plays a similar role as the increase of time. This is understood by considering that microorganisms with strong swimming ability move faster upstream and downstream, thus effectively requiring less time to achieve a similar longitudinal distribution for the vertical-mean concentration.
As shown in Figure 13, with the increase of time, the vertical concentration distribution C v becomes lower at the bottom and higher at the surface, reflecting a migration of the microorganisms to the water surface due to gravitaxis. However, with the increase of Pe s , the intensity of gravitactic accumulation decreases. Qualitatively, when a microorganism hit the free-water surface, it is reflected in an opposite direction and swims back to the water column, and then again, migrates to the water surface due to gravitaxis. If the swimming ability of microorganisms is strong, it takes a longer distance for them to reorient upwards after begin reflected by the water surface. Therefore, a larger Pe s can result in a thicker layer of the gravitactic accumulation in the vicinity of the water surface, and consequently a more uniform vertical concentration distribution. Note that the prescribed migration by Equation 8 is merely an idealization of the real-world complicated particle-wall interactions, as discussed in Section 2.2. It can only qualitatively characterize the influence of the presence of a surface boundary on the gravitactic accumulation. For the influence of the wall accumulation, we refer readers to the recent work of Wang et al. (2022b, Appendix B).

Influence of Shear Flow
Now we focus on the effect of flow shear on the dispersion of gyrotactic microorganisms. The strength of the background shear flow is quantified by Pe f , which we take three different values of 1, 5, and 10 in the analysis. The swimming Pèclet number Pe s and gyrotactic bias λ are fixed with values of 0.5 and 1, respectively. The same as discussed in previous sections, a patch of gyrotactic microorganisms is released instantaneously in the middle of the cross-section, and only the reflection boundary conditions are considered.

Drift Velocity and Dispersivity
From Equation 24 we see that −K 1 is the drift velocity. K 1 is negative for all shear strengths as illustrated in Figure 14a, describing the center of the microorganism cloud moving downstream. In addition, with a larger Pe f the value of K 1 is significantly greater, and the time needed to attain the steady state is also extended. This can be associated with the disruption of flow shear on the vertical migration of gyrotactic microorganisms. And the non-monotonic variations before approaching the long-time limit indicate the complex effect of boundary-microorganisms interaction.
The variation of K 2 in Figure 14b shows that the dispersion of the microorganism cloud is greatly enhanced by increasing Pe f . And similar to what is observed in Figure 14a for the K 1 -t curves, we see that the time to attain steady dispersion is obviously extended for a stronger shear flow. This phenomenon is in contrast to the classical dispersion of passive particles (Li et al., 2017), in which case the characteristic time to attain steady dispersion does not vary with the flow strength.

Concentration Distribution
As shown in Figure 15, microorganisms move downstream and migrate toward the water surface after being released. Specifically, with a larger Pe f the microorganism cloud moves downstream at a longer distance and longitudinally spreads more efficiently at the same time. Interestingly, we find that a large part of the microorganisms are temporally trapped in the vicinity of the vertical release position of z = 0.5 with Pe f = 10. This is due to gyrotactic trapping effect occurring in a region of high shear. However, the cloud eventually migrates to the water surface, because the gyrotactic trapping effect can only weaken instead of inhibiting the upward migration of the microorganisms.
A major difference in the longitudinal distribution of vertical-mean concentration, compared to that presented in previous sections, is that the symmetrical center of the Gaussian distribution no longer remains at the origin of x = 0, as can be observed in Figure 16. This is due to the convection and the temporal variation of the drift velocity of the cloud, which does not equal to the mean flow velocity during the transient transport process. When Pe f is larger, the convection of the background flow is stronger, and the center of the Gaussian distribution moves faster.
For the vertical concentration distribution as illustrated in Figure 17a, we see that as the increase of the strength of the shear flow, a high concentration area emerges in the middle part of the channel flow (Pe f = 10 at t = 5). Note that the sedimentary velocity of microorganisms is neglected in the current model. Thus, this gyrotactic trapping of microorganisms will gradually disappear and cannot be observed at a later time, for example, at t = 100, which is the same as we have discussed above. More detailed discussions can be found in the recent work of Wang et al. (2022b), who systematically studied the important role of sedimentation proposed by Ghorai (2016) in the form of thin trapping layers.

Effect of Turbulence on Microorganisms' Dispersion
The above theoretical analysis based on the Gill's generalized dispersion model is for the laminar open channel flow. That said, these observed characteristics of the microorganisms during transport in this idealized flow can qualitatively reflect their dispersion behavior in more complicated flow environment involving turbulence, since the fundamental elements for the dispersion remain the same which are either related to their own cell properties like the gyrotaxis and the swimming ability, or the vertical shear of the flow at a larger scale comparable to the characteristic length of the cross-section of the channel. On the other hand, it is also interesting to investigate quantitatively on how microorganisms transport in a more realistic flow can be different from that revealed by the above analysis. To be more realistic, a channel is better to be described in 3-dimensional, and turbulence is a necessary ingredient which needs to be included in the analysis.
Many previous studies (Durham et al., 2013;Zhan et al., 2014) have experimentally or numerically studied active transport in isotropic turbulent flows. Here, we primitively discuss the turbulent effect on microorganisms' dispersion by combining the direct numerical simulation (DNS) for the turbulent flow field and the random walk simulation for tracking the motions of the microorganisms, as that have been adopted previously for the dispersion of gyrotactic microorganisms  in vertical downwelling turbulent flows (Croze et al., 2013;Zhang et al., 2022). The advantage of directly resolving the turbulence by DNS is that we do not need to additionally adopt a turbulent diffusion model for the scattering of the transported microorganisms, the validation and calibration of which require a lot of experimental data that we currently do not have. Although the computational expense can be a challenge for high Reynolds number flows, here in this paper we discuss the contribution of turbulence on the transport of microorganisms as the turbulent intensity increases, which is affordable when the Reynolds number is not too high.

Direct Numerical Simulation
The fluid considered here is assumed to be incompressible and Newtonian, thus the fluid flow is governed by the following unsteady Navier-Stokes equations: where u designates the fluid velocity vector field. p f is the pressure. ρ f denotes the fluid density and ν f represents the kinematic viscosity.
As for the numerical methodology, the second-order Adams-Bashforth fractional-step scheme is used for the time advancement. All the terms are treated explicitly except for the vertical viscous term, for which the Crank-Nicolson scheme is utilized. The central finite-difference is utilized for spatial discretization with fourth-order accuracy in periodic directions and second-order accuracy in wall-normal direction. A staggered mesh configuration is employed for the strong coupling of velocity and pressure fields. More detailed descriptions about the numerical implementation can be found in our previous work (Gong & Fu, 2021;Gong et al., 2022Gong et al., , 2023.
The DNS of open-channel flow is performed in a rectangular domain of L x × L y × L z = 6πh × 2πh × h, where h is the channel height, and L x , L y , and L z denote the computational size in stream-wise, span-wise and wall-normal directions, respectively. Periodic conditions are imposed for both the carrier fluid and active cells in stream-wise and spanwise directions. No-slip condition is used for the fluid flow in the bottom wall and free-slip in the free-water-surface. In addition, completely elastic collisions are imposed when the distances between the cells and walls are smaller than the radius of the cell. The computational domain is discretized by 756, 360, 100 grids in stream-wise, spanwise and wall-normal direction, respectively. The corresponding dimensionless spatial resolutions for the turbulent cases are Δx*u τ /ν f = 4.5, Δy*u τ /ν f = 3.1, and Δz*u τ /ν f ∈ {0.16, 2.8}, where u τ is the frictional velocity. The time step is specified as Δt*D r = 2 × 10 −4 , and the corresponding CFL number is much smaller than unity.

Random Walk Simulation
Because we consider a dilute suspension, the cell-fluid interaction is neglected. Therefore, we can perform a similar random walk simulation (as described in Section 4) based on the velocity field at each moment while the fluid velocity is not affected by the presence of cells. Following Thorn and Bearon (2010), we use vector p to represent the swimming direction of a cell in a three-dimensional channel flow. In the simulation, the time is discrete and the corresponding dimensionless governing equation (Thorn & Bearon, 2010, Equations 25 and 26) for each time step δt can be written as where x is the position vector and U is the fluid velocity at the particle position, and U is calculated by the quadratic Lagrange interpolation method from the 27 closest fluid grids to the particle position. We have only considered the case of spherical cells and k is the unit vector of the vertical direction. √ 2 represents the white noise for the orientation vector p. Taking p as a polar axis, the resulted p′ disturbed by the noise follows a distribution, where the polar angle (of p′ with respect to p) θ′ ∼ N(0, 2δ t ) and the azimuthal angle ′ ∼  [0,2 ] . N and  represent the normal distribution and the uniform distribution, respectively.

Comparisons Between Dispersion in Laminar and Turbulent Flows
It is well known that in the turbulent flow, the small scale chaotic flow will eventually contribute to the scattering of the substances, in a way similar to that by a diffusion process as characterized by the turbulent diffusivity.
Here, to highlight the pure effect of the scattering by the chaotic flow, we simulate the two cases with different Reynolds numbers (Re m ∈ {286, 2860}) but with the same flow Péclet number Pe f . Namely, the mean flow rates are kept the same, while Re m = U m H*/ν* are controlled by changing the kinematic viscosity ν*. Note that Re m = 286 represents the laminar flow case, and Re m = 2860 is for the turbulent flow case.
At the very beginning, as a baseline demonstrating the predominance of the active behavior of the microorganisms, we adopt Pe f = 0.2 and Pe s = 0.5 for the simulation. We first present the velocity profiles in Figure 18 to clearly show the differences between the laminar and turbulent flows. In spite of these differences, we observe in Figure 19 that the two curves of vertical concentration distribution respectively for laminar and turbulent flows overlap. We note that as a demonstration for the limiting case of the transport process dominated by the swimming ability of the microorganisms, this result of identical vertical concentration distributions is expected since the contribution of the flow profile to the concentration distribution is negligible. That is, this specific distribution as observed in Figure 19 should be the same as the case even with no flow at all. More specifically, the dominance of swimming can also be quantitatively illustrated by the ratio of V s /U m = Pe s / Pe f = 2.5, meaning that the microorganisms can swim by a velocity that is more than twice the velocity of the mean flow, which definitely corresponds to a very weak flow or an extremely fast-swimming microorganism that may not exist. Another curve in Figure 19 is the result of analytical solution as obtained in previous sections, the comparison of it with the rest of the two curves in fact reveals the effect of channel geometry on the results: the 2-D case as theoretically investigated above, compared with the 3-D case considered in this section.
Then, we continue to consider the case with Pe f = 2 and Pe s = 0.5, in which case the contribution of the microorganisms' swimming plays a relatively important role during transport. By Figure 20, two curves respectively for laminar and turbulent flows no longer overlap but start to deviate from each other. The reduction of the surface accumulation for the turbulent flow case can be seen clearly by this vertical concentration distribution. Note that the angular velocity of the swimming direction of a microorganism is affected not only by the gravitational torque but also by the vorticity of the flows, as described by Equation 36. Thus, the upward swimming direction of cells will be more easily disturbed under such a chaotic flow by the random-like shear-induced rotation. Additionally, the turbulence-induced translational diffusion can also to some extent smooth the vertical concentration gradient, although this effect is not strong in this considered case.
A more detailed illustration regarding the concentration distributions in laminar and turbulent flows is shown in Figure 21. It is observed that the chaotic turbulent flow can suppress the surface accumulation and the dispersion of the gyrotactic microorganisms. Other than that, the mean velocity for the downstream motions of the microorganisms' cloud is slower, and the longitudinal scattering of it is seen to be weakened for the turbulent flow case. As a quantitative demonstration, we present the drift velocity and Taylor dispersivity for both the laminar and turbulent flows in Figure 22. The difference in the drift velocity is mainly due to the change of the velocity profile. For the Taylor dispersivity, note that the swimming-induced diffusion effect is dominant in the transport process because the advective effect is not strong. The chaotic turbulent flow can in fact enhance the effective rotational diffusion, and thus can further weaken the swimming-induced diffusion effect, resulting in a smaller K 2 .
The evolution of the longitudinal distribution of the cross-sectional mean concentration in Figure 23 shows the difference between the dispersion processes in laminar and turbulent flows. The concentration curve of the turbulent flow moves downstream slower and the peak concentration is higher in value. This phenomenon is characterized by the drift and dispersivity as shown in Figure 22 (see also the discussion therein). Generally, for the turbulent flow, since there is not as much proportion of microorganisms stayed close to the surface where they move downstream with the flow at the highest flow velocity, the overall drift velocity of the concentration cloud is slower. And a greater vertical mixing of the microorganisms corresponds to a  Finally, we consider the case with a much larger Pe f such that the contribution of the turbulence far more exceeds that of the active motions of the microorganisms (Pe f = 20, Pe s = 0.5). We focus on the vertical concentration and mean concentration, as shown in Figure 24. Although the turbulent disturbance has weakened the swimming-induced diffusion, the translational diffusion is greatly enhanced. The turbulent diffusivity is so large that the accumulation of the microorganisms at the vicinity of the free-water-surface is further weakened, with the overall vertical concentration distribution more close to a uniform distribution, although a greater concentration at the surface can still be observed. This is equivalent to active particles with a very large swimming ability which reduces the gyrotaxis effect, in which case the active particles act like the passive substances. With this large vertical diffusion effect, the necessary time for the initial non-uniform concentration distribution to approach the asymptotic uniform distribution also greatly decreased.

Concluding Remarks
In the present work, Gill's generalized dispersion model for passive substances is extended to investigate the transport of active particles. The effects of swimming, gyrotaxis, and shear strength on the dispersion of gyrotactic microorganisms in a laminar open-channel flow are discussed in detail. Through series expansion based on the longitudinal derivatives of cross-sectional mean concentration, we calculate the transient dispersion characteristics, specifically the drift velocity and the effective dispersion coefficient, based on which the evolution of the concentration distribution can be described.
Microorganisms with larger gyrotactic bias λ have a stronger tendency to form gravitactic accumulation at the vicinity of the free-water surface. For gyrotactic microorganisms, a non-uniform vertical distribution of the concentration can be observed, leading to great differences in microorganism concentration near the water surface and that at the bottom of the channel. In contrast, the vertical concentration distribution of non-gyrotactic microorganisms (λ = 0) eventually becomes uniform regardless of the initial condition of the transport. In addition, gravitactic accumulation reduces the longitudinal dispersion of microorganisms in the presence of no ambient fluid flows, since a larger gyrotactic bias means that the swimming of microorganisms is more projected to the direction of vertical instead of the longitudinal.
When microorganisms have a stronger swimming ability, they spread faster in all directions, resulting in larger longitudinal dispersivity K 2 and flatter longitudinal distribution of the vertical-mean concentration, under the circumstance of no ambient fluid flows. Another effect of the stronger swimming ability is that microorganisms take a longer distance to reorient upwards after being reflected by the boundary of the free-water surface, leading to a more uniform vertical concentration distribution.
With a stronger strength of the shear flow, the drift velocity and dispersion coefficient of microorganisms are significantly enhanced, and the time to approach the asymptotic dispersion regime is prolonged. Meanwhile, Figure 22. The evolution of (a) drift velocity K 1 and (b) dispersivity K 2 for laminar and turbulent flows. In all cases, λ = 1, Pe s = 0.5 and Pe f = 2.0. microorganisms move downstream faster and spread more efficiently in the longitudinal direction due to a stronger strength of advection and shear dispersion, respectively. And the center of mass, or the symmetrical center of the Gaussian distribution of the microorganisms can drift away from the origin of a moving frame with the mean velocity of the flow, demonstrating the transient initial effect of the transport when the drift velocity of the concentration cloud is not a constant. In a region of high shear, a large number of microorganisms can be temporarily trapped but eventually migrate to the water surface, which is understood as the gyrotactic trapping effect, which has been systematically investigated by previous studies (Durham et al., 2009;Ghorai, 2016;Wang et al., 2022b).
In this study, we have primarily considered the dispersion of gyrotactic microorganisms in laminar open-channel flows. As a demonstration on how turbulence will affect the observed transport characteristics as compared with that of the laminar flow, in Section 6 we have combined the techniques of DNS to resolve the turbulence and the random walk simulation to track the trajectories of microorganisms. We have adopted two Reynolds numbers to respectively describe the case of laminar (Re m = 286) and turbulent (Re m = 2,860) flows. Then, by considering different combinations of Pe f and Pe s which represent a variation of the relative importance of microorganisms' active swimming compared to the turbulence, we are able to reveal that as the increase of turbulence, the vertical concentration distribution will be more uniform, and the contribution of the active behavior of the microorganisms to the dispersion characteristics like the drift and Taylor dispersivity will decrease accordingly.
Some additional issues are worth exploring in the future as follows. First, we mainly consider the idealized reflective boundary conditions, which may not be able to fully capture the complex behaviors of microorganisms. It is necessary to take into account the interactions between microorganisms and the boundaries (either the free-water surface or the channel bottom), including the steric and hydrodynamic effects (Jakuszeit et al., 2019;Lushi et al., 2017). In addition, the interactions of microorganisms with the solid wall and the water surface may be very different, therefore, the influence of interface characteristics needs to be considered.
Second, we have only considered the dispersion of dilute microorganisms suspension. In practice, microorganisms impose stress on the fluid, which cannot be ignored in a relatively dense suspension (Saintillan, 2018). In order to study the dispersion of high concentration suspension, the interactions between microorganisms and fluid, and that among microorganisms should be considered.
Third, in this study, microorganisms are simplified as spheres, but most microorganisms are of the shape of prolate spheroid in practice. The increase of the aspect ratio of the single cell can affect the clustering and the average settling speed of the microorganisms by affecting the rotation rate and drag anisotropy of the microorganisms (Ardekani et al., 2017).
Finally, although good agreement is found between results obtained through Gill's generalized model and that by the random walk numerical simulation at relatively large times, the ability of Gill's generalized dispersion model in predicting the concentration distribution is not entirely satisfactory at the transient initial stage of the transport. This is due to limited information considered by the second-order generalized model (i.e., only the exact first-and second-order moments), excluding the non-Gaussian properties of the concentration cloud. Higher-order Gill's generalized model may be considered in further work.

Data Availability Statement
All the information for the theoretical analysis can be found in the text, and no experimental data are used. The code for performing the DNS in this paper is freely available at https://github.com/GongZheng-Justin/ CP3d.