# Should hydrodynamics be taken into account when calculating the growth rate of microalgae in a photobioreactor ?

J. Ignacio Fierro U., Liu-Di Lu, Olivier Bernard

<figure><img src="image_1.png" /></figure>

## To cite this version:

J. Ignacio Fierro U., Liu-Di Lu, Olivier Bernard. Should hydrodynamics be taken into account when calculating the growth rate of microalgae in a photobioreactor?. SIAM Journal on Applied Mathematics, 2025, 85 (4), pp.1-21. 10.1137/24M1693374. hal-05130895v2

HAL Id: hal-05130895

https://inria.hal.science/hal-05130895v2

Submitted on 27 Jun 2025

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# SHOULD HYDRODYNAMICS BE TAKEN INTO ACCOUNT WHEN CALCULATING THE GROWTH RATE OF MICROALGAE IN A PHOTOBIOREACTOR ?

J. IGNACIO FIERRO U. \*, LIU-DI LU \dagger OLIVIER BERNARD \ddagger

**Abstract.** Microalgae, as photosynthetic organisms, are cultivated in photobioreactors for various industrial applications. Light intensity, a critical factor influencing their growth rate, is inherently non-uniform within photobioreactors. In regions distant from the illuminated surface, microalgae experience photolimitation due to insufficient photon availability, hindering optimal activation of the photosynthetic machinery. Conversely, near the illuminated surface, excessive light intensity can damage key photosynthetic proteins, leading to photoinhibition. While mixing in photobioreactors does not alter the light gradient, it influences the light exposure history of cells through hydrodynamic advection. In this study, we employ Han's mechanistic model to describe the dynamics of photon harvesting and its consequences, including photoinhibition and photolimitation. First, we calculate the time-averaged growth rate for arbitrary continuous light signals, revealing how mixing impacts growth under the assumption of periodic light signals generated by hydrodynamics. Next, we address the computational challenge of estimating growth rates in photobioreactors using computational fluid dynamics (CFD), modeling a single-phase incompressible fluid. Finally, we analyze the case of a raceway pond, evaluating errors arising when growth rate is estimated without accounting for hydrodynamics. We demonstrate that the gain in growth is related to the cell movement along the light gradient. Our results show that in predominantly laminar hydrodynamic regimes, hydrodynamics has only a marginal effect on microalgal growth. Moreover, we show that the average productivity can be estimated based on a static approximation of the average growth rate taking into account the light distribution, with an error lower than 10%.

**Key words.** microalgae, mixing, computational fluid dynamics, raceway pond, Han model.

**MSC codes.** 92B99, 34A05, 34C25, 34A34, 76B07

**1. Introduction.** Microalgae are capable of converting CO₂ into biomass using energy from visible light. These microorganisms are cultivated industrially in either open or closed photobioreactors [19]. Open systems such as raceway ponds are simple and cost-effective, shallow, oval-shaped channels, mixed with a paddle wheel. Closed systems, like tubular reactors, represent more advanced technology with more controlled culturing conditions. These reactors can operate in batch mode or with continuous addition of growth medium. Two key characteristics define these systems: (1) microalgae cells act as light-absorbing particles, creating a heterogeneous light distribution within the reactor. Areas near the light source experience high illumination, while deeper regions remain in darkness. (2) Intensive mixing is employed to prevent biomass sedimentation and ensure uniform nutrient distribution. As a result, cells are advected through light gradients, experiencing alternating periods of high and low light intensity.

Photon harvesting in microalgae is a dynamic process, and the average growth rate in a photobioreactor emerges from the complex interplay between photosystem dynamics and hydrodynamics [7]. Accurately modeling this interaction is challenging, as it requires accurate representation of 1/ reactor hydrodynamics, 2/ light distribution within the reactor, 3/ dynamic response of photon harvesting in response to light variations. This approach, which explicitly accounts for the light history of cells, is

<footnote>\*UNIV. GRENOBLE ALPES, INRIA, CNRS, GRENOBLE INP, GIPSA-LAB, F-38000 GRENOBLE, FRANCE</footnote>

<footnote>†SECTION DE MATHÉMATIQUES, UNIVERSITÉ DE GENÈVE</footnote>

<footnote>‡CENTRE INRIA D'UNIVERSITÉ CÔTE D'AZUR, GREENOWL TEAM</footnote>
classified as type III in [3]. However, the complexity of type III models often limits
their practical application, leading to the use of simplified growth rate ($\mu$) calcula-
tions. Type I models depend solely on incident light at the reactor surface, while
type II models incorporate simple light transfer models like the Lambert–Beer law,
coupled with Monod-like or Haldane-type functions for $\mu(I)$. Unlike type I and II
models, which are static, type III models are dynamic.

Photosystem dynamics are typically described by three-population models, such
as the Han model [14] or the Eilers–Peeters model [11]. These models use three dif-
ferential equations to represent the probability of the photosystem being in one of
three states: open and ready to process photons, closed and processing photons,
or damaged due to excess light energy. These models capture two key phenom-
ena—photoabsorption and photoinhibition—which operate on different timescales.
The recovery rate after damage is significantly slower than photon harvesting, en-
abling a slow-fast approximation in the three-population model [15].

Cells in photobioreactor are exposed to continuous, fluctuating light signals over
time resulting from their trajectories in the light gradient. In practice, they are
most of the time exposed to suboptimal light conditions, reducing photosynthetic
efficiency. Recent studies have explored optimizing light absorption in these reactors,
by introducing specific topographies [5] or enhancing vertical mixing [6].

Using the Han model, we investigate the impact of fluctuating light signals on
growth rates. By considering the typical timescales of light variations, we simplify the
Han model using a slow-fast approximation and compare its predictions with those
assuming steady-state photosystems. We first apply our coupled model to simple
periodic light signals, commonly used in laboratory-scale photobioreactors. Then, we
reconstruct the cell trajectories in a raceway pond using computational fluid dynamics
(CFD) simulations. The trajectories of the cells in the light gradient provide realistic
light patterns whose effect on photosynthesis can be assessed using Han's model.

This paper is organized as follows. We first present the hydrodynamic and biolog-
ical models in Section 2, where we compare two strategies for computing the growth
rate: $\mu_A$ is a more realistic computation accounting for the dynamics of the photo-
systems, and $\mu_S$ is an approximation easier to compute assuming a static response of
the photosynthetic apparatus. The main results of this paper, is a characterization
of the relationship between the dynamic growth rate ($\mu_A$) and the static approxima-
tion ($\mu_S$). To better understand this relationship, we provide theoretical analysis in
Section 3 and study the impact of continuous periodic light signals. We demonstrate
that all solutions of the biological model converge to a unique periodic solution. Nu-
merical studies are provided in Section 4, where we first illustrate the mixing in a
photobioreactor using two typical periodic light signals. We then analyze the race-
way pond using a CFD model to simulate the motion of particles, tracking the light
perceived by individual microalgae. We compare the actual average growth rate and
the static approximation by taking space into account. A detailed discussion is given
in Section 5, where we comment our results and their applications on the design of
photobioreactors. Finally, we conclude by demonstrating that the average growth rate
computed assuming the steady state of the photosystems is a reliable approximation.

## 2. Hydrodynamic and biological models.

2.1. Computational fluid dynamic model and cell tracking. The water flow in a raceway pond can be simulated with CFD, which integrates the Navier–Stokes equations. Several studies have used CFD to simulate the velocity field in open ponds [25, 27]. Lagrangian approaches have also been used to assess the mixing effi-
<figure><img src="image_3.png" /><figcaption>Fig. 2.1: 2D mesh of the simulated raceway pond.</figcaption></figure>

ciency in algae cultures. For example, the mixing length is computed in [1] as a result of different paddle wheel velocities. In our study, we consider a real raceway pond from the Environmental Biotechnology Laboratory of INRAE Narbonne in France [17]. We use a layer-averaged Euler and Navier-Stokes model for the numerical simulation of incompressible free surface, as presented in [2]. The meshing of the raceway pond consists of a fixed 2D triangular mesh of the bottom of the raceway, as shown in Figure 2.1. The layers are defined by the water depth, giving the third dimension of the system. This approximation of the Navier-Stokes equations is more accurate than the well-known shallow water system. It employs a multi-layered model based on a Galerkin-type approximation of the velocity field, utilizing piecewise constant basis functions. A decisive advantage for raceway ponds is also that this discretization of Navier-Stokes equations represents the free surface (and the waves) more simply, yet more accurately than the classical approaches. The incompressible and hydrostatic Navier-Stokes system with free surface is given by

$$ (2.1) \qquad \begin{aligned} & \nabla \cdot \mathbf{U} = 0, \\ & \frac{\partial \mathbf{u}}{\partial t} + \nabla_{x,y} \cdot (\mathbf{u} \otimes \mathbf{u}) + \frac{\partial \mathbf{u}w}{\partial z} = \frac{1}{\rho_0} \nabla_{x,y} \cdot \sigma + \frac{\mu}{\rho_0} \frac{\partial^2 \mathbf{u}}{\partial z^2} + \mathbf{F}, \\ & \frac{\partial p}{\partial z} = -\rho_0 g, \end{aligned} $$

where $\mathbf{U} = (u, v, w)^T$ is the velocity of the liquid, $\mathbf{u} = (u, v)^T$ is the horizontal velocity, $\sigma$ is $-pI_d + \Sigma$, where $\Sigma = \mu \nabla_{x,y} \mathbf{u}$ is the total stress tensor, $p$ is the pressure, $g$ is the gravity acceleration constant, $\rho_0$ is the fluid density and $\mu$ is the viscosity coefficient. The fluid is assumed to be Newtonian. The hydrostatic Navier-Stokes system is relevant here, since vertical acceleration is negligible compared to the horizontal acceleration. The paddle wheel is indirectly represented by the force $\mathbf{F}$ as used in [4],

$$ (2.2) \qquad \mathbf{F} = F \left( \sqrt{(x - x_{\text{wheel}})^2 + (z - z_{\text{wheel}})^2} \omega \right)^2 \begin{pmatrix} \cos(\theta) \\ 0 \\ \sin(\theta) \end{pmatrix}, $$

where $F$ is a constant, $\theta$ is the angle between the vertical axis and the blade, $\omega = \dot{\theta}$, and $x_{\text{wheel}}$, $z_{\text{wheel}}$ are the coordinates of the paddle wheel in the $x$ and $z$ axis respectively. Note that the force does not affect the $y$-axis, which is parallel to the central axis of
the paddle wheel. It has been shown that using a 2D (horizontal and vertical axis) representation of the raceway hydrodynamics [4] is computationally more efficient than using model (2.1)-(2.2).

The fluid domain along the z-axis is delimited by the free surface denoted by $η(t, x, y)$, then the system (2.1) is completed with the following boundary condition:

$$ (2.3) \qquad \frac{\partial \eta}{\partial t} + \mathbf{u} \cdot \nabla_{x,y} \eta = 0. $$

The CFD model was validated in [2] using a velocity sensor in a smaller raceway pond. Our simulations are very similar to the one in [17], carried out with a commercial software, where a physical model of the paddle wheel was implemented.

We assume that the microalgae have the same density as the medium, so that the trajectories of the cells match that of the background flow. Under this hypothesis, Lagrangian's trajectories of these particles can be reconstructed from the Eulerian description (2.1) denoted by $(\mathbf{X}_n)_{i=1}^{N_{\text{par}}}$, where $N_{\text{par}}$ is the number of simulated particles. The position of each particle $\mathbf{X}_n(t) = (x_n, y_n, z_n)^T$ is computed by solving the equation

$$ (2.4) \qquad \frac{d\mathbf{X}_n(t)}{dt} = \mathbf{U}(t), \quad \mathbf{X}_n(0) = (x_{n0}, y_{n0}, z_{n0})^T, $$

where $(x_{n0}, y_{n0}, z_{n0})$ is the initial position of the particle. We denote by $\Omega$ the domain of the raceway. We initiate the position of the particles randomly, following a uniform distribution in the domain $\Omega$.

Similar works on the same pond [16, 23] have considered a growth model tracking the position of Lagrangian trajectories, using the same model (2.1), (2.2), and (2.4). In our study, we have considerably improved the numerical schemes, so that much more particle trajectories could be simulated with a higher accuracy. In particular, it is challenging to reproduce with the simulation the expected equidistribution of the cells along time, and much of the biases appeared in previous schemes were reduced. To illustrate this photobioreactor, we show a simulated raceway pond in Figure 2.2 together with the distribution of the velocity magnitude. The paddle wheel is positioned just above the red surface. Some cells (in white) are represented together with their respective streamlines. The trajectories of the cells are mainly horizontal in the straight sections of the raceway pond.

**2.2. Light distribution within the raceway pond.** We assume that a population of microalgae cells, each with a radius of 10 (μm), follows fluid streamlines within the light gradient. Due to the incompressibility of the fluid and the cells' density being equal to that of the medium, the particles are, in theory, uniformly distributed throughout the raceway. This assumption justifies the use of the Lambert–Beer law to approximate the light distribution, taking into account the scattering and absorption effects of the cells. The vertical motion of the cells in the light gradient can therefore provide the individual light history of each cell. We assume that light strikes perpendicularly to the ground surface. To reach a cell $\mathbf{X}_n$, the light must travel a distance $η(t, x_n, y_n) - z_n$, from the free surface. Then the light signal $I_n$ perceived by this particle is computed with the Lambert–Beer law as

$$ (2.5) \qquad I_n(t) = I_0 e^{-\xi(\eta(t, x_n, y_n) - z_n)}, $$

where $I_0$ is the light perceived at the free surface, $\xi > 0$ is the light extinction constant.
<figure><img src="image_4.png" /><figcaption>Fig. 2.2: Representation of a simulated raceway pond and the velocity field of some particles. The height corresponds to 0.3 (m). **Left:** Simulated particles representation (paddle wheel simulated at 12.5 RPM). **Top Right:** Lagrangian tracers of few particles (paddle wheel simulated at 20 RPM). **Bottom Right:** The streamline of some particles’ trajectories (paddle wheel simulated at 15 RPM).</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>Fig. 2.3: **Left:** Illustration of the evolution between the three states in the Han model. **Right:** Parameters of the Han model used in this study.</figcaption></figure>

<table><thead><tr><td>Parameter</td><td>Value [Unit]</td></tr></thead><tbody><tr><td>k<sub>r</sub></td><td>6.8 · 10<sup>-3</sup> [s<sup>-1</sup>]</td></tr><tr><td>k<sub>d</sub></td><td>2.99 · 10<sup>-4</sup> [-]</td></tr><tr><td>τ</td><td>0.25 [s]</td></tr><tr><td>σ<sub>P</sub></td><td>0.047 [m<sup>2</sup>μmol<sup>-1</sup>]</td></tr><tr><td>α</td><td>8.7 · 10<sup>-6</sup> [-]</td></tr></tbody></table>

**2.3. Biological model.** To describe how photons are harvested by the photosystems under a varying light intensity, we consider the mechanistic model of Han [14], which characterizes the process of photon harvesting with possible photoinhibition induced by the photodamage of the photosystem II (PSII). There are three possible states for PSII: open or reactive state A, closed or activated state B, and inhibited or damaged state C. The relation of these three states are schematically presented in Figure 2.3. The dynamics of PSII can be described by the differential equations:

$$ (2.6) \quad \frac{dA}{dt} = -I\sigma_{P}A + \frac{B}{\tau}, \quad \frac{dB}{dt} = I\sigma_{P}A - \frac{B}{\tau} + k_{r}C - k_{d}\sigma_{P}IB, \quad \frac{dC}{dt} = -k_{r}C + k_{d}\sigma_{P}IB, $$

where $\sigma_P$ (μmol$^{-1}$m$^2$) is the effective cross-section of PSII, $I$ (μmolm$^{-2}$s$^{-1}$) is the light intensity perceived by the microalgae, $\tau$ (s) is the minimal time required for an electron to transfer from water on the donor side of the photosynthetic unit to the terminal electron acceptors, which is also called the turnover time, $k_d$ (–) is the damage rate and $k_r$ (s$^{-1}$) is the recovery rate of PSII. A, B and C represent the probability distribution of each state, therefore

$$ (2.7) \qquad A + B + C = 1. $$
The algal growth rate is assumed to be proportional to the open state A and the light intensity I. More precisely, the growth rate given by the kinetic model (2.6) corresponds to

$$ (2.8) \qquad \mu_A(I, A) := \alpha \sigma_P I A, $$

where $\alpha$ is a constant of proportionality relating the flux of electrons triggered by the photons and cell growth. At steady state, the state A of the system (2.6) tends towards $A_S$ (see, e.g., [7]). This corresponds exactly to a Haldane model where the growth rate is a function of light, that is

$$ (2.9) \qquad \mu_S(I) = \alpha \sigma_P I A_S = \alpha \sigma_P I \frac{1}{1 + \tau \sigma_P I + \frac{k_d}{k_r} \tau (\sigma_P I)^2}. $$

The maximum of $\mu_S$ is given by $\mu_{max} = \frac{\alpha}{\tau + 2\sqrt{\frac{k_d}{k_r}}\tau}$, and it is reached when the light intensity value is

$$ (2.10) \qquad I_{\text{opt}} = \frac{1}{\sigma_P \sqrt{\frac{k_d}{k_r} \tau}}. $$

As shown in [12], using (2.7) and substituting it into the system (2.6), we end up with two equations

$$ (2.11) \qquad \begin{aligned} \frac{dA}{dt} &= -\left(\sigma_P I + \frac{1}{\tau}\right)A - \frac{1}{\tau}C + \frac{1}{\tau}, \\ \frac{dC}{dt} &= k_d \left[-\sigma_P I A - \left(\sigma_P I + \frac{k_r}{k_d}\right)C + \sigma_P I\right]. \end{aligned} $$

In practice, the factor $k_d$ is in the range of $10^{-4}$, as shown in Figure 2.3. Therefore, system (2.11) has slow/fast timescales, and we can consider the slow manifold proposed in [18, Chapter 11], where A rapidly reaches a pseudo steady state depending on the value of C:

$$ (2.12) \qquad A = \frac{1-C}{1+\tau\sigma_P I}, $$

This reduces the system (2.11) to a single equation of C,

$$ (2.13) \qquad \frac{dC}{dt} = -(\gamma(I) + k_r)C + \gamma(I), $$

with $\gamma(I) := \frac{k_d \tau (\sigma_P I)^2}{1+\tau \sigma_P I}$. In particular, the steady state of C is given by $C_S := \frac{\gamma(I)}{\gamma(I)+k_r} = \frac{\frac{k_d}{k_r} \tau (\sigma_P I)^2}{1+\tau \sigma_P I + \frac{k_d}{k_r} \tau (\sigma_P I)^2}$. Consequently, the steady state of A can be obtained by substituting $C_S$ into (2.12).

For a continuous bounded light signal $I: [0, +\infty) \to [I_{min}, I_{max}]$ with $0 \le I_{min} \le I_{max}$, we present an estimation of the actual growth rate $\mu_A$ as a function of the static approximation $\mu_S$ assuming the Han model in equilibrium. To simplify the notation, we write $\mu_A(t)$ instead of $\mu_A(I(t), A(I(0); t))$ and $\mu_S(t)$ instead of $\mu_S(I(t))$. We now demonstrate that the initial condition rapidly does not affect the dynamics anymore and can therefore be neglected after a time interval of $10/k_r$ (about 25 minutes), much smaller than the typical timescale of growth.
**PROPOSITION 2.1.** Let $I: [0, +\infty) \to [I_{\min}, I_{\max}]$ be a continuous bounded light signal. Assuming that $A(0) = 0$, the growth rate $\mu_A$ can be written as

$$ (2.14) \qquad \mu_A(t) = \mu_S(t)(\gamma(I(t)) + k_r) \int_0^t e^{-\int_s^t \gamma(I(w)) + k_r \,dw} \,ds. $$

Furthermore, the following estimations hold

$$ (2.15) \qquad \mu_S(t) \frac{\gamma(I(t)) + k_r}{\gamma(I_{\max}) + k_r} \left(1 - e^{-(\gamma(I_{\max}) + k_r)t}\right) \le \mu_A(t) \le \mu_S(t) \frac{\gamma(I(t)) + k_r}{\gamma(I_{\min}) + k_r}. $$

In general, for any initial condition $A(0) = A_0$, the associated growth rate converges to (2.14) when $t$ goes to infinity.

*Proof.* The general solution of (2.13) is given by

$$ (2.16) \qquad C(t) = C(0)e^{-\int_0^t \gamma(I(w)) + k_r \,dw} + \int_0^t \gamma(I(s))e^{-\int_s^t \gamma(I(w)) + k_r \,dw} \,ds. $$

Using then (2.12), the state A can be written as

$$ (2.17) \qquad \begin{aligned} A(A(0); t) = &A(0) \frac{1 + \tau \sigma_P I(0)}{1 + \tau \sigma_P I(t)} e^{-\int_0^t \gamma(I(w)) + k_r \,dw} \\ &+ A_S(t)(\gamma(I(t)) + k_r) \int_0^t e^{-\int_s^t \gamma(I(w)) + k_r \,dw} \,ds, \end{aligned} $$

which depends on the initial value $A(0)$. If $A(0) = 0$, the latter becomes,

$$ (2.18) \qquad A(0; t) = A_S(t)(\gamma(I(t)) + k_r) \int_0^t e^{-\int_s^t \gamma(I(w)) + k_r \,dw} \,ds. $$

Then, equation (2.14) is deduced from the definitions (2.8) and (2.9). The upper bound in (2.15) follows from the fact that $\gamma$ is an increasing function of $I$ and

$$ \begin{aligned} \int_0^t e^{-\int_s^t \gamma(I(w)) + k_r \,dw} \,ds &\le \int_0^t e^{-\int_s^t \gamma(I_{\min}) + k_r \,dw} \,ds = \frac{1 - e^{-(\gamma(I_{\min}) + k_r)t}}{\gamma(I_{\min}) + k_r} \\ &\le \frac{1}{\gamma(I_{\min}) + k_r}. \end{aligned} $$

For a similar reason, we find $\int_0^t e^{-\int_s^t \gamma(I(w)) + k_r \,dw} \,ds \ge \int_0^t e^{-\int_s^t \gamma(I_{\max}) + k_r \,dw} \,ds = \frac{1 - e^{-(\gamma(I_{\max}) + k_r)t}}{\gamma(I_{\max}) + k_r}$. Substituting these two inequalities into (2.18) and using the definitions (2.8) and (2.9), we obtain the estimation (2.15). Finally, if $A(0) = A_0 \ne 0$, using (2.17) and (2.18), we have

$$ (2.19) \qquad \begin{aligned} |A(A(0); t) - A(0; t)| &= |A(0)| \frac{1 + \tau \sigma_P I(0)}{1 + \tau \sigma_P I(t)} e^{-\int_0^t \gamma(I(w)) + k_r \,dw} \\ &\le |A(0)| \frac{1 + \tau \sigma_P I_{\max}}{1 + \tau \sigma_P I_{\min}} e^{-k_r t}. \end{aligned} $$

The latter converges to zero at a rate $k_r$, thus $A(A(0); t)$ converges to $A(0; t)$. $\square$

Based on Proposition 2.1, for large timescales, we will assume that the initial condition of the state A is zero.
**2.4. Growth rate in the raceway pond.** To estimate the growth rate in the raceway pond, we need to define the average growth rate of all simulated particles moving within the photobioreactor. We first define the time-averaged dynamic growth rate and its static approximation by

$$ (2.20) \qquad \bar{\mu}_A := \frac{1}{T} \int_0^T \mu_A(t) \, dt, \quad \bar{\mu}_S := \frac{1}{T} \int_0^T \mu_S(t) \, dt. $$

Then the space-time-averaged dynamic growth rate and its static approximation are defined as

$$ (2.21) \qquad \bar{\bar{\mu}}_A = \frac{1}{N_{\text{par}}} \sum_{n=1}^{N_{\text{par}}} \bar{\mu}_A(I_n), \qquad \bar{\bar{\mu}}_S = \frac{1}{N_{\text{par}}} \sum_{n=1}^{N_{\text{par}}} \bar{\mu}_S(I_n). $$

In the Eulerian description, the static approximation of the growth rate $\mu_S$ can be computed for each point in the raceway domain $\Omega$. At each point $(x, y, z) \in \Omega$, the perceived light is computed using (2.5) as $I(x, y, z) = I_0 e^{-\xi(\eta(t,x,y)-z)}$, and the volume-averaged growth rate can be defined by

$$ (2.22) \qquad \bar{\bar{\mu}}_{\Omega} = \frac{1}{V(\Omega)} \int_{\Omega} \mu_S(I(x, y, z)) \, dxdydz, $$

where $V(\Omega)$ is the volume of the raceway. In [16], the static approximation $\bar{\bar{\mu}}_S$ is used to compute the growth rate in a raceway pond. In [29], the relationship between $\bar{\bar{\mu}}_S$ and $\bar{\bar{\mu}}_{\Omega}$ is discussed in a different type of photobioreactor. We focus on a more accurate computation of the growth rate using $\bar{\bar{\mu}}_A$ and characterize the gain compared to the static approximation computed with $\bar{\bar{\mu}}_S$.

Instead of analyzing the flow of microalgae cells at a fixed location (Eulerian perspective), tracking Lagrangian trajectories provides a natural approach to monitor the light exposure perceived by each cell. The light history derived from these trajectories is then used to compute the growth rate $\bar{\bar{\mu}}_A$. Notably, the growth rate $\bar{\bar{\mu}}_A$ cannot be defined within a Eulerian framework, as the concept of light history is lost in this perspective.

**3. Formal analysis of the average growth rate.**

**3.1. Time-averaged growth rate.** We first assess the difference between the average growth rate $\bar{\bar{\mu}}_A$ and the approximation $\bar{\bar{\mu}}_S$ defined in (2.20).

**THEOREM 3.1.** Let $I$ be a bounded continuous light signal perceived by a single cell. Assuming that $A(0) = 0$, for a given time period $T$, the time-averaged growth rate $\bar{\bar{\mu}}_A$ can be written as the sum of the time-averaged static approximation $\bar{\bar{\mu}}_S$, and a correction term $\bar{\bar{\mu}}_H$ representing the dynamical gain due to mixing:

$$ (3.1) \qquad \bar{\bar{\mu}}_A = \bar{\bar{\mu}}_S + \bar{\bar{\mu}}_H + \mathcal{O}(1/T), $$

where

$$ (3.2) \qquad \bar{\bar{\mu}}_H := \frac{1}{T} \int_0^T \frac{d\mu_S(t)}{dt} \phi(t) \, dt, \quad \phi(t) := \int_0^t e^{-\int_s^t \gamma(I(w))+k_r \, dw} \, ds, $$

and $\mathcal{O}(1/T) = -\frac{\mu_S(T)}{T} \phi(T)$ which goes to zero when $T \to +\infty$.
**Interpretation of Theorem 3.1:** In a system mainly laminar, where cells stay at a constant depth (or are slowly advected vertically), for which the growth rate gain due to vertical mixing $\bar{\mu}_H$ is negligible, the average growth rate over a sufficiently long time period $T$ can be accurately computed using the static approximation $\bar{\mu}_S$ based on the Haldane model. In a system with high velocities in the direction of the light gradient, this approximation must be refined.

*Proof.* Using (2.14) and integration by parts, we find

$$ 
\begin{aligned}
\int_0^T \mu_A(t) dt &= \int_0^T \mu_S(t)(\gamma(I(t)) + k_r) \int_0^t e^{-\int_s^t \gamma(I(w)) + k_r \, dw} \, ds \, dt \\
&= \int_0^T \int_s^T \mu_S(t)(\gamma(I(t)) + k_r) e^{-\int_s^t \gamma(I(w)) + k_r \, dw} \, dt \, ds \\
&= \int_0^T \mu_S(s) - \mu_S(T) e^{-\int_s^T \gamma(I(w)) + k_r \, dw} \, ds \\
&\quad + \int_0^T \int_s^T \frac{d\mu_S(t)}{dt} e^{-\int_s^t \gamma(I(w)) + k_r \, dw} \, dt \, ds.
\end{aligned}
$$

Dividing the latter by $T$, we get (3.1). Moreover, for a given time period $T$, the function $\phi$ defined in (3.2) can be upper bounded by

$$ (3.3) \qquad \phi(T) \le \int_0^T e^{-(\gamma(I_{\min})+k_r)(T-s)} \, ds = \frac{1 - e^{-(\gamma(I_{\min})+k_r)T}}{\gamma(I_{\min}) + k_r}, $$

and lower bounded by

$$ (3.4) \qquad \phi(T) \ge \int_0^T e^{-(\gamma(I_{\max})+k_r)(T-s)} \, ds = \frac{1 - e^{-(\gamma(I_{\max})+k_r)T}}{\gamma(I_{\max}) + k_r}. $$

This reveals the fact that $\mathcal{O}(1/T) \to 0$ as $T \to \infty$. $\square$

<figure><img src="image_6.png" /></figure>

Note that $|\mathcal{O}(1/T)| \le \frac{\mu_{\max}}{T} \frac{1}{\gamma(I_{\min})+k_r} \le \frac{\mu_{\max}}{T} \frac{1}{k_r}$. Then, for $T \gg 1/k_r$, we can approximate $\bar{\mu}_A$ by $\bar{\mu}_S + \bar{\mu}_H$, where $\bar{\mu}_H$ is defined in (3.2). The quantity $1/k_r$ corresponds to the time needed to go from the state $C$ to $B$ in the Han model, and $T \gg 1/k_r$ means that the period $T$ must be large enough to incorporate the effect of recovery.

**3.2. Space-time-averaged growth rate.** We discuss here how the hydrodynamics of the raceway pond affects the growth rate $\bar{\mu}_A$ defined in (2.21). For each particle, the velocity in the $z$-axis is given by (2.4), i.e., $v_z(t, \mathbf{X}_n) = \mathbf{U}_w(\mathbf{X}_n(t))$. The next proposition relates the velocity of the $z$-axis and the growth rate $\bar{\mu}_A$.

**PROPOSITION 3.2.** *The space-time-averaged dynamic growth rate can be bounded by the average velocity $v_z$ in the $z$-axis as*

$$ (3.5) \qquad \bar{\mu}_A \le \bar{\mu}_S + \frac{\alpha \xi \sigma_P I_0}{k_r N_{\text{par}}} \sum_{n=1}^{N_{\text{par}}} \frac{1}{T} \int_0^T |v_z(t, \mathbf{X}_n)| \, dt + \frac{\mu_{\max}}{Tk_r}. $$

*Proof.* For each particle $\mathbf{X}_n$ and its perceived light signal $I_n$, we have $|\bar{\mu}_H(I_n)| \le \frac{1}{T} \int_0^T \left| \frac{d\mu_S(I_n(t))}{dt} \right| \phi(I_n(t)) \, dt = \frac{1}{T} \int_0^T \left| \frac{d\mu_S(I_n(t))}{dI_n} \frac{dI_n(t)}{dt} \right| \phi(I_n(t)) \, dt$. Using then the in-
equality (3.3), we obtain

$$ \begin{aligned} |\bar{\mu}_H(I_n)| &\le \frac{1}{T} \int_0^T \left| \frac{d\mu_S(I_n(t))}{dI_n} \frac{dI_n(t)}{dt} \right| \frac{(1 - e^{-(\gamma(I_{n,\min}+k_r)t})}{\gamma(I_{n,\min}) + k_r} dt, \\ &\le \frac{1}{Tk_r} \int_0^T \left| \frac{d\mu_S(I_n(t))}{dI_n} \frac{dI_n(t)}{dt} \right| dt, \end{aligned} $$

where $I_{n,\min}$ is the minimum value of $I_n$ in the interval $[0, T]$. From the Lambert-Berth law, we have $\frac{d}{dt}I_n(t) = -\xi I_0 e^{-\xi(\eta(t,x_n,y_n)-z_n)} (\frac{\partial\eta}{\partial t} + \frac{\partial\eta}{\partial x}v_x + \frac{\partial\eta}{\partial y}v_y - v_z)$. Using (2.3) in the above, the derivative w.r.t. the light is $\frac{d}{dt}I_n(t) = -\xi I_n v_z(t, \mathbf{X}_n)$. It follows that $|\bar{\mu}_H(I_n)| \le \frac{\xi}{Tk_r} \int_0^T \left| \frac{d\mu_S(I_n(t))}{dI_n} I_n(t)v_z(t, \mathbf{X}_n) \right| dt$. Note that $\frac{d}{dI}\mu_S(I) = \frac{\alpha\sigma_P(1-\frac{k_d}{k_r}\tau(\sigma I)^2)}{(1+\tau\sigma_P I+\frac{k_d}{k_r}\tau(\sigma I)^2)^2} \le \alpha\sigma_P$. Using the latter and the fact that $I_n(t) \le I_0$, we have $|\bar{\mu}_H(I_n)| \le \frac{\alpha\xi\sigma_P I_0}{Tk_r} \int_0^T |v_z(t, \mathbf{X}_n)| dt$. Using (3.1), we get

$$ \begin{aligned} \bar{\mu}_A &= \bar{\mu}_S + \frac{1}{N_{\text{par}}} \sum_{n=1}^{N_{\text{par}}} \bar{\mu}_H(I_n) - \frac{\mu_S(I_n(t))}{T} \phi(I_n(t)), \\ &\le \bar{\mu}_S + \frac{1}{N_{\text{par}}} \sum_{n=1}^{N_{\text{par}}} |\bar{\mu}_H(I_n)| + \frac{\mu_S(I_n(t))}{T} \phi(I_n(t)), \\ &\le \bar{\mu}_S + \frac{1}{N_{\text{par}}} \sum_{n=1}^{N_{\text{par}}} \frac{\alpha\xi\sigma_P I_0}{Tk_r} \int_0^T |v_z(t, \mathbf{X}_n)| dt + \frac{\mu_{\text{max}}}{Tk_r}, \end{aligned} $$

which proves the inequality (3.5). □

**Interpretation of Proposition 3.2:** We can ignore the term $\mu_{\max}/T$ in (3.5) when the time period $T$ is large enough. Then, the difference between the space-time-averaged growth rate $\bar{\mu}_A$ and $\bar{\mu}_S$ cannot be greater than $\frac{\alpha\xi\sigma_P I_0}{k_r}|v_z|$ with $|v_z| := \frac{1}{N_{\text{par}}} \sum_{n=1}^{N_{\text{par}}} \int_0^T |v_z(t, \mathbf{X}_n)| dt$. Here, $|v_z|$ is an indicator of how mixed the photobioreactor is on the z-axis. A non-mixed photobioreactor will have $|v_z| = 0$, leading to no difference between $\bar{\mu}_S$ and $\bar{\mu}_A$. Furthermore, the value of the light extinction $\xi$ also plays an important role. Lower values of this parameter present less variations of the light gradient inside the reactor.

**3.3. Periodic light signals.** When the light signal is periodic, we can be more accurate in the results presented in Theorem 3.1. Let $T$ now be the period of the continuous light signal $I$ perceived by the cell, we have the following result.

**PROPOSITION 3.3.** Let $I$ a continuous periodic function, i.e., $I(t + T) = I(t)$, $\forall t \in [0, +\infty)$. Then, all solutions of (2.13) converge to a unique periodic solution.

*Proof.* The periodic solution $C_p$ is obtained by imposing $C_p(T) = C_p(0)$ in (2.16), and we find $C_p(0) = \frac{\int_0^T \gamma(I(s))e^{-\int_s^T \gamma(I(w)) + k_r dw} ds}{1 - e^{-\int_0^T \gamma(I(w)) + k_r dw}}$. Let $C$ be a solution of (2.13) and $\varepsilon = C - C_p$, we have $\dot{\varepsilon} = -(\gamma(I) + k_r)\varepsilon$. The solution is given by $\varepsilon(t) = \varepsilon(0)e^{-\int_0^t \gamma(I(s)) ds}e^{-k_r t}$, which can be bounded by $|\varepsilon(t)| \le |\varepsilon(0)|e^{-(\gamma(I_{\min})+k_r)t}$. Then $\varepsilon$ goes to zero as $t$ goes to $\infty$. □

**PROPOSITION 3.4 (Dynamic growth rate for periodic signal).** Let us consider a continuous periodic light signal $I$ of period $T$, i.e., $I(T) = I(0)$. Then the dynamic
growth rate $\mu_A$ associated with the only periodic solution $C_p$ is given by:

$$ (3.6) \qquad \mu_A(t) = \mu_S(t)(\gamma(I(t)) + k_r) \left[ \phi(t) + \frac{\phi(T)O(T)}{1 - O(T)} \right]. $$

where $\phi$ is given by (3.2) and $O(T) = e^{-\int_0^T \gamma(I(w)) + k_r \, dw}$.

*Proof.* From (2.17), we have

$$ \mu_A(T) = \mu_S(T)(\gamma(I(t)) + k_r)\phi(T) + \mu_A(0)\frac{I(T)}{I(0)}\frac{1 + \tau\sigma_P I(0)}{1 + \tau\sigma_P I(T)}e^{-\int_0^T \gamma(I(w)) + k_r \, dw}. $$

As the function $I$ is periodic $I(T) = I(0)$, then $C_p(T) = C_p(0)$, and the growth rate is also periodic. Imposing $\mu_A(T) = \mu_A(0)$, we obtain $\mu_A(T) = \mu_S(T)(\gamma(I(t)) + k_r)\phi(T) + \mu_A(T)e^{-\int_0^T \gamma(I(w)) + k_r \, dw}$, and thus $\mu_A(T) = \mu_S(T)\frac{(\gamma(I(t)) + k_r)\phi(T)}{1 - e^{-\int_0^T \gamma(I(w)) + k_r \, dw}}$. Substituting then into (2.14), we find (3.6). $\square$

In the same way, we obtain the following result for the time-averaged dynamic growth rate $\bar{\mu}_A$ in the periodic case.

**THEOREM 3.5 (Time-averaged growth rate in the periodic case).** *Considering a periodic light signal $I$ of period $T$. The time-averaged dynamic growth rate $\bar{\mu}_A$ is given by:*

$$ (3.7) \qquad \bar{\mu}_A = \bar{\mu}_S + \frac{1}{1 - O(T)} \left[ \bar{\mu}_H + \frac{O(T)}{T} \int_0^T \phi(T) - \phi(t) \, dt \right]. $$

**Interpretation of Theorem 3.5:** This theorem clarifies the results of Theorem 3.1, and better characterizes the dynamical component of the growth rate $\bar{\mu}_H$ which must be added to the static approximation $\bar{\mu}_S$ when the velocities along the light gradient are marked.

*Proof.* The proof follows the same steps as in Theorem 3.1 using the dynamic growth rate computed in the periodic case (3.6). $\square$

Note that when $T$ is large enough, meaning that $T \gg 1/k_r$, $O(T) \to 0$, leading again to the same approximation in the non-periodic case, i.e., $\bar{\mu}_A \approx \bar{\mu}_S + \bar{\mu}_H$.

**3.4. Characterization of the dynamical component $\bar{\mu}_H$ as a function of the light regime.** We have already seen that the actual time-averaged growth rate can be approximated by $\bar{\mu}_A \approx \bar{\mu}_S + \bar{\mu}_H$, and $\bar{\mu}_H + O(1/T)$ goes to zero as the period $T$ goes to $\infty$. In other words, $\bar{\mu}_S$ becomes an accurate approximation for large period $T$, i.e., for slow movement against the light gradient. In this section, we provide a lower and upper bounds for $\bar{\mu}_H$ to understand its relationship with the variation rate of the light signal $I$.

**PROPOSITION 3.6.** Let $I: [0, T] \to [0, +\infty)$ be a light signal of class $C^1$, such that $I_{\min} \le I(t) \le I_{\max}$, $\forall t \in [0, T]$. Assume that all the stationary points of the function $t \mapsto \mu_S(I(t))$ are isolated. For $\{t_i\}_{i=0}^m$ a partition of $[0, T]$, where $t_0 = 0$, $t_m = T$, and $\frac{d\mu_S(I(t_i))}{dt} = 0$, then

$$ (3.8) \qquad -\frac{\mu_{\max}}{T} \frac{|\mathcal{I}_-|}{\gamma(I_{\max}) + k_r} \le \bar{\mu}_H \le \frac{\mu_{\max}}{T} \frac{|\mathcal{I}_+|}{\gamma(I_{\min}) + k_r}, $$
where $|\mathcal{I}_+|$ and $|\mathcal{I}_-|$ are the cardinality of the sets:

$$ \mathcal{I}_+ := \left\{ i : \frac{d\mu_S(I(t))}{dt} > 0 \ \forall t \in (t_i, t_{i+1}) \right\}, \quad \mathcal{I}_- := \left\{ i : \frac{d\mu_S(I(t))}{dt} < 0 \ \forall t \in (t_i, t_{i+1}) \right\}. $$

*Proof.* According to the definition of $|\mathcal{I}_+|$ and $|\mathcal{I}_-|$, $\bar{\mu}_H$ can be decomposed into positive and negative parts as $\bar{\mu}_H = \bar{\mu}_H^+ + \bar{\mu}_H^-$, where $\bar{\mu}_H^+ := \frac{1}{T} \sum_{i \in \mathcal{I}_+} \int_{t_i}^{t_{i+1}} \frac{d\mu_S(t)}{dt} \phi(t) dt$ and $\bar{\mu}_H^- := \frac{1}{T} \sum_{i \in \mathcal{I}_-} \int_{t_i}^{t_{i+1}} \frac{d\mu_S(t)}{dt} \phi(t) dt$. Now, we can give an upper bound for the positive part $\bar{\mu}_H^+$ using (3.3),

$$ \bar{\mu}_H^+ \le \frac{1}{T} \frac{1}{\gamma(I_{\min}) + k_r} \sum_{i \in \mathcal{I}_+} \int_{t_i}^{t_{i+1}} \frac{d\mu_S(t)}{dt} dt = \frac{1}{T} \frac{1}{\gamma(I_{\min}) + k_r} \sum_{i \in \mathcal{I}_+} \mu_S(t_{i+1}) - \mu_S(t_i). $$

In the same way, we can give a lower bound for the negative part using (3.4),

$$ \bar{\mu}_H \ge \frac{1}{T} \frac{1}{\gamma(I_{\max}) + k_r} \sum_{i \in \mathcal{I}_-} \int_{t_i}^{t_{i+1}} \frac{d\mu_S(t)}{dt} dt = \frac{1}{T} \frac{1}{\gamma(I_{\max}) + k_r} \sum_{i \in \mathcal{I}_-} \mu_S(t_{i+1}) - \mu_S(t_i). $$

Note that $\bar{\mu}_H^- \le \bar{\mu}_H \le \bar{\mu}_H^+$ due to the sign of each term. Then, using the lower bound of $\bar{\mu}_H^-$ and the upper bound of $\bar{\mu}_H^+$, we find $\frac{1}{T} \frac{\sum_{i \in \mathcal{I}_-} \mu_S(t_{i+1}) - \mu_S(t_i)}{\gamma(I_{\max}) + k_r} \le \bar{\mu}_H \le \frac{1}{T} \frac{\sum_{i \in \mathcal{I}_+} \mu_S(t_{i+1}) - \mu_S(t_i)}{\gamma(I_{\min}) + k_r}$. Then, as $\mu_S(t_{i+1}) - \mu_S(t_i) \le \mu_{\max}$, we have $-\frac{1}{T} \frac{\sum_{i \in \mathcal{I}_-} \mu_{\max}}{\gamma(I_{\max}) + k_r} \le \bar{\mu}_H \le \frac{1}{T} \frac{\sum_{i \in \mathcal{I}_+} \mu_{\max}}{\gamma(I_{\min}) + k_r}$. $\square$

**Interpretation of Proposition 3.6:** A particle in a constantly mixed photo-bioreactor cannot remain at the same depth. Therefore, the perceived light intensity $I$ cannot be constant in any interval, and the function $t \mapsto \mu_S(I(t))$ can only have isolated stationary points. The value of $\bar{\mu}_H$ is bounded by the number of times the derivative of $\mu_S(I(t))$ changes its sign. If $I$ is a periodic function as the one analyzed in Section 3.3, then the two sums $\sum_{i \in \mathcal{I}_+} \mu(t_{i+1}) - \mu(t_i)$ and $\sum_{i \in \mathcal{I}_-} \mu(t_{i+1}) - \mu(t_i)$ are independent of $T$. And if $T \to +\infty$, then $\bar{\mu}_H$ converges to 0.

An efficient trajectory for a microalgae maximizes $\bar{\mu}_H$, i.e., increases growth rate due to fast changes of light along the light gradient. So, the question is which type of perceived light signals provide a higher value of $\bar{\mu}_H$. To address this point, we give a lower bound for $\bar{\mu}_H^+$:

$$ \begin{aligned} \bar{\mu}_H^+ &= \frac{1}{T} \sum_{i \in \mathcal{I}_+} \int_{t_i}^{t_{i+1}} \frac{d\mu_S(I(t))}{dt} \int_0^t e^{-\int_s^t \gamma(I(w)) + k_r} dw \, ds dt \\ &\ge \frac{1}{T} \sum_{i \in \mathcal{I}_+} \int_{t_i}^{t_{i+1}} \frac{d\mu_S(I(t))}{dt} \int_0^{t_i} e^{-(\gamma(I_{\max}) + k_r)(t_{i+1} - s)} \, ds dt, \\ &= \frac{1}{T} \sum_{i \in \mathcal{I}_+} \frac{\delta_i}{\gamma(I_{\max}) + k_r} (\mu_S(I(t_{i+1})) - \mu_S(I(t_i))), \end{aligned} $$

where $\delta_i = e^{-(\gamma(I_{\max}) + k_r)(t_{i+1} - t_i)} (1 - e^{-(\gamma(I_{\max}) + k_r)t_i})$.

The quantity $\bar{\mu}_H$ is instrumental to understand how growth $\bar{\mu}_A$ is stimulated in a photobioreactor by the frequent oscillations along the light gradient. This key observation will be determinant for the design and operation of photobioreactors. If
we compute $\bar{\mu}_H$ for different mixing strategies, it will inform about the most efficient mixing strategy for growth. An increase in the value of $\bar{\mu}_H$ results from higher $\bar{\mu}_H^+$ and $\bar{\mu}_H^-$ . To increase the value of $\bar{\mu}_H^+$, we can look at the value of $\delta_i(\mu_S(I(t_{i+1}))-\mu_S(I(t_i)))$. The value of $\delta_i$ is larger for a shorter time interval ($t_i, t_{i+1}$), i.e., for faster movements along the light gradient. Then, each short interval where the value $\mu_S(I(t_i))$ moves to a higher value $\mu_S(I(t_{i+1}))$ helps to increase the value of $\bar{\mu}_H^+$, and consequently, increase the value of $\bar{\mu}_A$. Although computing $\bar{\mu}_H$ is challenging, the previous observation offers a key insight to optimize productivity in a photobioreactor. In the particular case of the raceway pond, Equation (3.8) can be used to get the bound:

$$
-\frac{\mu_{\max}}{T} \frac{1}{N_{\text{par}}} \sum_{n=1}^{N_{\text{par}}} \frac{|\mathcal{I}_{-}(I_n)|}{\gamma(I_{\max}) + k_r} \le \frac{1}{N_{\text{par}}} \sum_{n=1}^{N_{\text{par}}} \bar{\mu}_{H}(I_n) \le \frac{\mu_{\max}}{T} \frac{1}{N_{\text{par}}} \sum_{n=1}^{N_{\text{par}}} \frac{|\mathcal{I}_{+}(I_n)|}{\gamma(I_{\min}) + k_r}.
$$

The quantity $\frac{1}{N_{\text{par}}} \sum_{n=1}^{N_{\text{par}}} |\mathcal{I}_-(I_n)|$ represents the average number of visits below the depth $z_{\text{opt}}$. In the same way, the quantity $\frac{1}{N_{\text{par}}} \sum_{n=1}^{N_{\text{par}}} |\mathcal{I}_+(I_n)|$ represents the average number of potential visits above the optimal depth $z_{\text{opt}}$ (see e.g., Figure 4.4). The above bounds are not tight and they provide a potential for growth enhancement if vertical hydrodynamics is sufficiently efficient.

4. Numerical study.

**4.1. Periodic sketchy examples.** For the numerical test, we chose the pa-
rameter values of the Han model from [13] as shown in Figure 2.3. To illustrate the
behavior of the dynamic growth rate, and specifically to compare it with the static
one, we consider two examples of periodic light signals. We first consider a simple
periodic function for the depth of a cell

$$
(4.1) \qquad z(t) = \frac{H_0}{2} \left( 1 + \sin \left( \frac{2\pi}{T} t \right) \right).
$$

The light signal, using the Lambert–Beer law, is given by $I(t) = I_0 e^{-\xi z(t)}$ and the static approximation of the growth rate is $\mu_S(t) = \mu_S(I(t))$. In this case, the value of $\bar{\mu}_S$ is independent of $T$ using a change of variable $s = t/T$,

$$
\bar{\mu}_S = \frac{1}{T} \int_0^T \mu_S(I(t)) dt = \int_0^1 \frac{\alpha \sigma_P I_0 e^{-\xi z(s)}}{1 + \tau \sigma_P I_0 e^{-\xi z(s)} + \frac{k_d}{k_r} \tau (\sigma_P I_0 e^{-\xi z(s)})^2} ds.
$$

In Figure 4.1a, $\mu_A$ is illustrated for the periodic solution given by (3.6). As expected, for $T$ large enough, the difference $\bar{\mu}_A - \bar{\mu}_S$ can be approximated by $\bar{\mu}_H$. Numerically, the value of $\bar{\mu}_H$ is close to zero. In fact, their relative difference in percentage, computed as $100 \times \frac{\bar{\mu}_A - \bar{\mu}_S}{\bar{\mu}_S}$, is lower than 2%, as illustrated in Figure 4.1b. As shown in the same figure, $\bar{\mu}_A$ is always greater than $\bar{\mu}_S$ and the difference between them becomes smaller as $T$ increases. When light varies slower, the approximation $\bar{\mu}_S$ for $\bar{\mu}_A$ is still accurate.

The actual growth rate $\mu_A$ is not always greater than the static approximation $\mu_S$ as it is shown in the second example,

$$
(4.2) \qquad z(t) = H_0 - 4H_0 \frac{e^{-\left(\sin\left(\frac{2\pi t}{T}\right)\right) - \frac{1}{2}}}{\left(1 + e^{-\left(\sin\left(\frac{2\pi t}{T}\right)\right) - \frac{1}{2}}\right)^2},
$$
<figure><img src="image_7.png" /><figcaption>(a) The growth rate $\mu_A$ (in red) and the static (blue) time-averaged growth rate $\bar{\mu}_A$ (contin- approximation $\mu_S$ (in blue) for two different uous red line) and static approximation $\bar{\mu}_S$ values of $T$. As $T$ increases, $\mu_A$ approaches (blue line) for different values of $T$. The to $\mu_S$. dashed green line shows the difference in per- centage.</figcaption></figure>

Fig. 4.1: Dynamic growth rate and static approximation comparison when the light signal is given by (4.1).

which also has a period $T$. This function describes a movement closer to the surface and does not travel to the deepest part of the culture as in (4.1). By doing the same change of variable as above, the time-averaged static approximation $\bar{\mu}_S$ is still independent of $T$. Similarly to Figure 4.1b, Figure 4.2b also shows the convergence of $\bar{\mu}_A$ to $\bar{\mu}_S$ when $T$ goes to $\infty$. The difference in percentage between $\bar{\mu}_S$ and $\bar{\mu}_A$ is greater than in the first example. Moreover, unlike in Figure 4.1b, the value of $\bar{\mu}_A$ in Figure 4.2b is always lower than the value of $\bar{\mu}_S$. For both examples, we use $H_0 = 0.3\text{m}$ and $I_0 = 1200\text{m}^{-2}\mu\text{mol}$.

As predicted in Proposition 3.3, for both examples, the estimation of the actual growth rate $\mu_A$ converges to $\mu_S$ (see Figure 4.1a and Figure 4.2a) as $T \to \infty$. Mixing effects are seen when the period is lower than 1500 (s). The first example shows how mixing can increase the growth rate, while the second can decrease it. Figure 4.3 presents the cell trajectories in each case. In the second example, the trajectory remains above the optimal depth, i.e., in the photoinhibited part of the reactor. In this case, the effects of photoinhibition are enhanced by mixing, which explains the decrease in growth rate. However, the growth rate in the second example is still greater than in the first example. Above all, the second case represents a trajectory staying above the optimal light for growth, which means that, in the photobioreactor, other trajectories will always stay below this threshold, so that the resulting overall average growth rate is likely to be low. In conclusion, it is difficult to find an optimal mixing pattern that meets several criteria at the same time. In addition, due to the incompressibility of the fluid, leading to particle equidistribution, the trajectories of all cells fully explore the space. Only CFD simulation can eventually reconstruct realistic light patterns.
<figure><img src="image_8.png" /><figcaption>(a) The growth rate $\mu_A$ (in red) and the static (blue) time-averaged growth rate $\bar{\mu}_A$ (contin- approximation $\mu_S$ (in blue) for two different uous red line) and static approximation $\bar{\mu}_S$ values of $T$. As $T$ increases, $\mu_A$ approaches (blue line) for different values of $T$. The dashed green line shows the difference in per- centage.</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>Fig. 4.3: Cell trajectories for two typical examples. The continuous line represents the trajectory given by (4.1) and the dashed line represents the trajectory associated with (4.2).</figcaption></figure>

**4.2. CFD simulations in the raceway pond.** To account for more realistic light signals, we investigate light driven by the hydrodynamics. Figure 4.4 illustrates the tracking of a single cell within a simulated raceway pond. We use the *Freshkiss3D* Python library for the CFD simulation. We distinguish two areas within the photobioreactor: the photoinhibited section (light green) where the light perceived is greater than the optimal light $I_{opt}$ defined in (2.10), and the photolimited section (dark green) where the light perceived is lower than $I_{opt}$. We denote by $z_{opt}$ the optimal depth, which corresponds to the depth at which the algae perceive the light $I_{opt}$ given by $z_{opt} = \frac{1}{\xi} \ln(\frac{I_0}{I_{opt}})$. When the particle travels from the photolimited
<figure><img src="image_10.png" /><figcaption>Fig. 4.4: Actual growth rate $\mu_A$ and static approximation $\mu_S$ for a simulated trajectory in a raceway pond operated at 20 RPM with an initial water height 0.3 (m). The maximum value of the function $\mu_S$ is denoted by $\mu_{max}$ and the depth at which this value is reached by $z_{opt}$.</figcaption></figure>

area to the photoinhibited area, crossing $z_{opt}$, the difference between $\mu_S$ and $\mu_A$ is more observable. This occurs when the particle moves faster on the z-axis than the photoinhibition mechanism. In this case, the actual growth rate $\mu_A$ could reach larger values than the estimation from the static approximation $\mu_S$. When the cell does not go through this section fast enough, then $\mu_S$ can accurately approximate $\mu_A$.

To account for the entire raceway system, we simulate 8 different hydrodynamical conditions by changing the velocity of the paddle wheel for a simulated time of one hour. The initial positions of the particles ($x_{n0}, y_{n0}, z_{n0}$) are randomly generated following an independent uniform distribution $x_{n0} \sim U[x_{min}, x_{max}], y_{n0} \sim U[y_{min}, y_{max}], z_{n0} \sim U[z_{min}, z_{max}]$, where the volume is defined by $B := [x_{min}, x_{max}] \times [y_{min}, y_{max}] \times [z_{min}, z_{max}]$ such that $\Omega \subset B$, i.e., we ignore the particles outside the domain $\Omega$ of the raceway pond. Due to this process, the simulations have a slightly different number of simulated particles $N_{par}$. The table on the right of Figure 4.5 shows the number of simulated particles for each simulation and the total CPU time. The software *Freshkiss3D* solves (2.1) to get the velocity field, at the same time it tracks the trajectories of particles by solving (2.4). The experiments were carried out on a computer with an Intel Xeon w-2223 processor running at 1200 MHz with a total of 15677 MB of RAM and Fedora version 39.

The set of particles necessary to estimate accurately $\bar{\bar{\mu}}_A$ and $\bar{\bar{\mu}}_S$ should be representative of the entire raceway pond, meaning that the distribution should be uniform, at least on the z-axis. As experimentally validated and confirmed by simulations, perfect mixing is reached in the raceway, even at the slowest velocities. In theory, due to the incompressibility of the fluid, all biochemical quantities should be uniformly distributed. An initial uniform distribution is imposed here, but after some time, the numerical error tends to accumulate turn after turn, and particles' distribution becomes higher both in the upper and lower layers. The particles' distribution progressively drifts from equidistribution to a distribution with two pits. This error is
<figure><img src="image_11.png" /><figcaption>Fig. 4.5: **Left:** Growth rate of the raceway simulated for 8 different velocities of the paddle wheel in Revolutions Per Minutes (RPM). In red, the growth rate $\overline{\bar{\mu}}_A^w$ and in blue, its static approximation $\overline{\bar{\mu}}_S^w$. In black, the growth rate computed with the volume defined by the free surface $\overline{\bar{\mu}}_\Omega$. In green, the difference between them in percentage. **Right:** Computation time in hrs and the number of simulated particles for each velocity.</figcaption></figure>

<table><thead><tr><th>RPM</th><th>N<sub>par</sub></th><th>Computation time (hrs)</th></tr></thead><tbody><tr><td>10</td><td>3094</td><td>714.8</td></tr><tr><td>12.5</td><td>3057</td><td>620.4</td></tr><tr><td>15</td><td>3064</td><td>742.2</td></tr><tr><td>17.5</td><td>3073</td><td>663.1</td></tr><tr><td>20</td><td>3064</td><td>779.3</td></tr><tr><td>22.5</td><td>3034</td><td>695.1</td></tr><tr><td>25</td><td>3034</td><td>810.5</td></tr><tr><td>27.5</td><td>3048</td><td>730.4</td></tr></tbody></table>

difficult to track, it results from the numerical approximations for solving Navier-Stokes equations, together with the numerical error when tracking each single cell in the Lagrangian approach. More particles could be simulated to mitigate this bias, but this is very CPU intensive and time consuming as shown in Figure 4.5.

Another approach consists in correcting the particles' distribution, so that the set of particles stays equidistributed. We propose a weighted average for correcting the numerical bias in the particles' distribution and more accurately compute $\overline{\bar{\mu}}_A$ and $\overline{\bar{\mu}}_S$,

$$ \overline{\bar{\mu}}_A^w = \frac{1}{T} \int_0^T \frac{\sum_{n=1}^{N_{\text{par}}} \frac{\mu_A(I_n(t))}{h(t, \mathbf{X}_n)}}{\sum_{n=1}^{N_{\text{par}}} \frac{1}{h(t, \mathbf{X}_n)}} dt \quad \text{and} \quad \overline{\bar{\mu}}_S^w = \frac{1}{T} \int_0^T \frac{\sum_{n=1}^{N_{\text{par}}} \frac{\mu_S(I_n(t))}{h(t, \mathbf{X}_n)}}{\sum_{n=1}^{N_{\text{par}}} \frac{1}{h(t, \mathbf{X}_n)}} dt. $$

Here, $\mu_A(I_n(\cdot))$ and $\mu_S(I_n(\cdot))$ are computed from the light history $I_n(\cdot)$ perceived by the algae, using (2.14) and (2.9) respectively, and $h(t, \cdot)$ is the probability density function of the depth of the particles at the instant $t$. If the free surface is perfectly flat with depth $H_0$, we have that $h(0, \cdot) = 1/H_0$ by construction (since all particles are initialized with uniform distribution on the z-axis), and we recover the expression (2.21). Figure 4.5 shows the computation of $\overline{\bar{\mu}}_A^w$ and $\overline{\bar{\mu}}_S^w$.

**4.3. Approximation of the growth rate.** When the particle number $N_{\text{par}}$ is large enough, the space-time-averaged static approximation defined in (2.21) converges to the volume-averaged growth rate (2.22) in the sense that

$$ (4.3) \qquad \lim_{N_{\text{par}} \to \infty} \overline{\bar{\mu}}_S(I_n) = \lim_{N_{\text{par}} \to \infty} \overline{\bar{\mu}}_S^w(I_n) = \overline{\bar{\mu}}_\Omega, $$

when the fluid is incompressible [29]. Figure 4.5 shows $\overline{\bar{\mu}}_\Omega$ computed for eight different paddle wheel velocities. We observe that this quantity is a constant independent of the simulated velocities. In fact, the volume-averaged growth rate $\overline{\bar{\mu}}_\Omega$ depends only
on the aquatic volume in the photobioreactor. Since the volume is constant in our numerical tests, this quantity can then be computed only considering the geometry of the raceway pond and the height of the water, making $\overline{\mu}_{\Omega}$ an efficient approximation of the growth rate. Two raceway ponds with equivalent shape lead to the same volume-averaged growth rate $\overline{\mu}_{\Omega}$ [20, Theorem 3]. The same holds for the corrected space-time-averaged static approximation $\overline{\mu}_{S}^{w}$ when $N_{par} \to \infty$. And in the limit case, $\overline{\mu}_{S}^{w}$ no longer sees the effect of hydrodynamics. This explains the reason why $\overline{\mu}_{S}^{w}$ has a small variation when the paddle wheel velocity changes in Figure 4.5. More precisely, the values of the corrected space-time-averaged static approximation $\overline{\mu}_{S}^{w}$ are 1.1395 (d⁻¹) for 10 RPM and 1.12080 (d⁻¹) for 27.5 PRM. The difference between the actual growth rate computed using the dynamic description $\overline{\mu}_{A}^{w}$ and the static approximation $\overline{\mu}_{S}^{w}$ is lower than 8% (at 27.5 RPM). Then, the actual growth rate $\overline{\mu}_{A}^{w}$ is almost not sensitive to the simulated velocity of the paddle wheel. In this way, it is not worthwhile to simulate hydrodynamics to obtain a more accurate measurement of the growth rate, as CFD simulations are very time consuming (see the table in Figure 4.5). Hydrodynamics have a minor effect on the overall growth rate, this conclusion is firmly related to the system presented here. A different photobioreactor could lead to different conclusions.

**5. Discussion.** Simulating the coupling between hydrodynamics and photosynthesis is a significant scientific challenge, involving multiple timescales and nonlinear dynamic models. In this study, we demonstrated how the average growth rate of microalgae in a bioreactor can be optimized through careful management of hydrodynamics. A key observation is that the increase in productivity is directly linked to the movement of cells along the light gradient. This dynamic movement enhances the growth rate compared to a static scenario, where cells remain immobile.

Better understanding hydrodynamics and their favorable impact on microalgae productivity paves the way for optimizing photobioreactors. In particular, it enables the design of systems capable of generating targeted movements along the light gradient, a key factor for maximizing microalgae growth. Other types of movements, while essential for ensuring mixing homogeneity, have no significant effect on productivity. Therefore, a strategic approach aimed at exploiting beneficial movements within the light gradient could significantly enhance the efficiency of microalgae cultivation. However, modeling this coupling between physics and biology presents several challenges. On one hand, the dynamics of photosynthesis in response to light fluctuations are complex to capture. Current models, which focus on the efficiency of photosystems, require more robust experimental validation [10]. On the other hand, hydrodynamic models have their own limitations, particularly when adopting a Lagrangian approach. Indeed, the property of iso-distribution of particles tends to degrade over time, an aspect often overlooked in previous studies but crucial in the context of periodic fluid circulations. To address this issue, we propose a strategy to compensate distribution biases, essential to avoid numerical drifts and erroneous conclusions.

Our study focuses on a simplified case, where we assume no cell sedimentation (i.e., cell density equal to that of the fluid) and perfect adherence to streamlines. A range of biological mechanisms such as photoacclimation and photoprotection with specific pigments to cope with high light have been neglected [10]. Additionally, we neglect external light fluctuations, such as those induced by day-night cycles. Although our approach is general from a mathematical perspective, the simulations are specifically applied to high-rate ponds. Similar studies for photobioreactors could validate these mechanisms in faster hydrodynamic regimes [26], where the productivity
gain from cell agitation could be significantly higher than that observed in raceway
ponds. Mixing can also indirectly promote growth for other reasons. On top of
avoiding sedimentation, it enhances the mass transfer with gases [8], leading to higher
CO2 transfer rate or increasing O2 outgassing, both being favorable for photosynthesis.
Here, the study was carried out assuming a constant medium turbidity, i.e., a constant
microalgal biomass concentration. A higher growth rate would probably support
a higher biomass in the reactor, reducing the average light in the reactor. This
effect should be further studied with dedicated models also taking into account other
timescale of the photosynthesis [24].

An important conclusion of this study is that the productivity gain associated
with mixing velocity stays in the range of a few percent of the static approximation
of the growth rate. This most probably explains why, despite recurrent statements
in the literature [28], the increase of productivity with more intense agitation has
never been clearly demonstrated experimentally. Given the high computational cost
of simulating reactor hydrodynamics and tracking particles to compute the average
growth rate, a correction factor applied to the growth parameter, accounting for the
hydrodynamics, is probably the most efficient numerical approach. Such correction
term should be calibrated on real systems, to automatically capture the photosynthesis
stimulation due to the cell movement in the light field. It is crucial to balance the
gain in productivity with the energy required for mixing, which increases with the
cube of the fluid velocity [9, 30]. Excessive agitation could lead to disproportionate
energy costs without a significant improvement in productivity [21]. These results
could serve as a basis for a life cycle assessment (LCA) to determine the optimal
agitation intensity, taking into account environmental and energy impacts [22]. Such
an approach would help identify a trade-off between productivity and sustainability.

**6. Conclusions.** We presented a numerical method for calculating the average growth rate in a photobioreactor, incorporating the effects of hydrodynamics. Our approach was based on the Han model to represent the photosynthesis dynamics and account for the light history of the microalgae. By coupling the Han model with hydrodynamics, we captured the influence of mixing devices on microalgae growth. We demonstrated analytically that cell advection along the light gradient is the determining reason to stimulate photosynthesis efficiency.

We simulated the light harvesting model within the light field generated by the
hydrodynamics of a raceway pond. Many publications have been dedicated to improve
the mixing in photobioreactors, but the criterion of the resulting algal productivity
was not taken into account. Our study opens new routes to more directly optimize
bioreactor productivity through hydrodynamics management. The relationship be-
tween the growth rate in the raceway and the vertical velocity (z-axis) can be extended
to other photobioreactors, depending on their hydrodynamics and internal light distri-
bution. This insight provides a basis for optimizing the reactor geometry and mixing
to maximize the benefits of hydrodynamic effects by ensuring cell movements along
the light gradient.

REFERENCES

[1] ALI, H., CHEEMA, T. A., YOON, H.-S., DO, Y., AND PARK, C. W. Numerical prediction of algae cell mixing feature in raceway ponds using particle tracing methods. *Biotechnology and bioengineering* **112**, 2 (2015), 297–307.

[2] ALLGEYER, S., BRISTEAU, M.-O., FROGER, D., HAMOUDA, R., JAUZEIN, V., MANGENEY, A., SAINTE-MARIE, J., SOUILLÉ, F., AND VALLÉE, M. Numerical approximation of the 3D
hydrostatic Navier-Stokes system with free surface. *ESAIM: Mathematical Modelling and Numerical Analysis* 53, 6 (2019), 1981–2024.

[3] BÉCHET, Q., SHILTON, A., AND GUIEYSSE, B. Modeling the effects of light and temperature on algae growth: state of the art and critical assessment for productivity prediction during outdoor cultivation. *Biotechnology advances* 31, 8 (2013), 1648–1663.

[4] BERNARD, O., BOULANGER, A.-C., BRISTEAU, M.-O., AND SAINTE-MARIE, J. A 2D model for hydrodynamics and biology coupling applied to algae growth simulations. *ESAIM: Mathematical Modelling and Numerical Analysis* 47, 5 (2013), 1387–1412.

[5] BERNARD, O., LU, L.-D., SAINTE-MARIE, J., AND SALOMON, J. Topography optimization for enhancing microalgal growth in raceway ponds. Submitted, 2024.

[6] BERNARD, O., LU, L.-D., AND SALOMON, J. Optimal periodic resource allocation in reactive dynamical systems: Application to microalgal production. *International Journal of Robust and Nonlinear Control* 33, 9 (2023), 4989–5010.

[7] BERNARD, O., MAIRET, F., AND CHACHUAT, B. Modelling of microalgae culture systems with applications to control and optimization. *Microalgae Biotechnology* (2015), 59–87.

[8] CAIA, M., BERNARD, O., AND BÉCHET, Q. Optimizing CO₂ transfer in algal open ponds. *Algal research* 35 (2018), 530–538.

[9] CHIARAMONTI, D., PRUSSI, M., CASINI, D., TREDICI, M. R., RODOLFI, L., BASSI, N., ZITTELLI, G. C., AND BONDIOLI, P. Review of energy balance in raceway ponds for microalgae cultivation: Re-thinking a traditional system is possible. *Applied Energy* 102 (2013), 101–111.

[10] DEMORY, D., COMBE, C., HARTMANN, P., TALEC, A., PRUVOST, E., HAMOUDA, R., SOUILLÉ, F., LAMARE, P.-O., BRISTEAU, M.-O., SAINTE-MARIE, J., ET AL. How do microalgae perceive light in a high-rate pond? towards more realistic lagrangian experiments. *Royal Society open science* 5, 5 (2018), 180523.

[11] EILERS, P., AND PEETERS, J. A model for the relationship between light intensity and the rate of photosynthesis in phytoplankton. *Ecological modelling* 42, 3-4 (1988), 199–215.

[12] FIERRO ULOA, J. I., LU, L.-D., AND BERNARD, O. Theoretical growth rate of microalgae under high/low-flashing light. *Journal of Mathematical Biology* 86, 4 (2023), 48.

[13] GRENIER, J., LOPES, F., BONNEFOND, H., AND BERNARD, O. Worldwide perspectives of rotating algal biofilm up-scaling. Submitted paper, 2020.

[14] HAN, B.-P. A mechanistic model of algal photoinhibition induced by photodamage to photosystem-II. *Journal of theoretical biology* 214, 4 (2002), 519–527.

[15] HARTMANN, P., BÉCHET, Q., AND BERNARD, O. The effect of photosynthesis time scales on microalgae productivity. *Bioprocess and biosystems engineering* 37, 1 (2014), 17–25.

[16] HARTMANN, P., DEMORY, D., COMBE, C., HAMOUDA, R., BOULANGER, A.-C., BRISTEAU, M.-O., SAINTE-MARIE, J., SIALVE, B., STEYER, J.-P., RABOUILLE, S., ET AL. Growth rate estimation of algae in raceway ponds: A novel approach. *IFAC Proceedings Volumes* 47, 3 (2014), 6216–6221.

[17] HREIZ, R., SIALVE, B., MORCHAIN, J., ESCUDIÉ, R., STEYER, J.-P., AND GUIRAUD, P. Experimental and numerical investigation of hydrodynamics in raceway reactors used for algaculture. *Chemical Engineering Journal* 250 (2014), 230–239.

[18] KHALIL, H. *Nonlinear Systems*. Pearson Education. Prentice Hall, 2002.

[19] MAGALHÃES, I. B., FERREIRA, J., DE SIQUEIRA CASTRO, J., DE ASSIS, L. R., AND CALIJURI, M. L. Agro-industrial wastewater-grown microalgae: A techno-environmental assessment of open and closed systems. *Science of The Total Environment* 834 (2022), 155282.

[20] MARTÍNEZ, C., MAIRET, F., AND BERNARD, O. Theory of turbid microalgae cultures. *Journal of theoretical biology* 456 (2018), 190–200.

[21] MENDOZA, J., GRANADOS, M., DE GODOS, I., ACIÉN, F., MOLINA, E., BANKS, C., AND HEAVEN, S. Fluid-dynamic characterization of real-scale raceway reactors for microalgae production. *Biomass and Bioenergy* 54 (2013), 267–275.

[22] MORALES, M., COLLET, P., LARDON, L., HÉLIAS, A., STEYER, J.-P., AND BERNARD, O. Life-cycle assessment of microalgal-based biofuel. *Biofuels from algae* (2019), 507–550.

[23] NIKOLAOU, A., BOOTH, P., GORDON, F., YANG, J., MATAR, O., AND CHACHUAT, B. Multi-physics modeling of light-limited microalgae growth in raceway ponds. *IFAC-PapersOnLine* 49, 26 (2016), 324–329.

[24] NIKOLAOU, A., HARTMANN, P., SCIANDRA, A., CHACHUAT, B., AND BERNARD, O. Dynamic coupling of photoacclimation and photoinhibition in a model of microalgae growth. *Journal of theoretical biology* 390 (2016), 61–72.

[25] PANDEY, R., AND PREMALATHA, M. Design and analysis of flow velocity distribution inside a raceway pond using computational fluid dynamics. *Bioprocess and Biosystems engineering* 40, 3 (2017), 439–450.
[26] PERNER, I., POSTEN, C., AND BRONESKE, J. CFD optimization of a plate photobioreactor used for cultivation of microalgae. *Engineering in life sciences* 3, 7 (2003), 287–291.

[27] PRUSSI, M., BUFFI, M., CASINI, D., CHIARAMONTI, D., MARTELLI, F., CARNEVALE, M., TREDICI, M. R., AND RODOLFI, L. Experimental and numerical investigations of mixing in raceway ponds for algae cultivation. *Biomass and bioenergy* 67 (2014), 390–400.

[28] RICHMOND, A., Ed. *Handbook of microalgal culture: biotechnology and applied phycology*, vol. 577. Wiley Online Library, 2004.

[29] TONG, Z.-X., LI, M.-J., YAN, J.-J., AND GU, Z.-L. A theoretical analysis of the hydrodynamic influence on the growth of microalgae in the photobioreactors with simple growth kinetics. *International Journal of Heat and Mass Transfer* 158 (2020), 119986.

[30] ZHU, C., JI, Y., DU, X., KONG, F., CHI, Z., AND ZHAO, Y. A smart and precise mixing strategy for efficient and cost-effective microalgae production in open ponds. *Science of the Total Environment* 852 (2022), 158515.