# Optimal bacterial resource allocation strategies in batch processing

Agustín Gabriel Yabo, Jean-Baptiste Caillau, Jean-Luc Gouzé

► To cite this version:

Agustín Gabriel Yabo, Jean-Baptiste Caillau, Jean-Luc Gouzé. Optimal bacterial resource allocation strategies in batch processing. SIAM Journal on Applied Mathematics, 2024, pp.S567-S591. 10.1137/22M1506328. hal-03710681v2

HAL Id: hal-03710681

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

Submitted on 25 Sep 2023

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished 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_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# OPTIMAL BACTERIAL RESOURCE ALLOCATION STRATEGIES IN BATCH PROCESSING*

AGUSTÍN GABRIEL YABO†, JEAN-BAPTISTE CAILLAU‡, AND JEAN-LUC GOUZÉ§

**Abstract.** The study of living microorganisms using resource allocation models has been key in elucidating natural behaviors of bacteria, by allowing allocation of microbial resources to be represented through optimal control strategies. The approach can also be applied to research in microbial cell factories, to investigate the optimal production of value-added compounds regulated by an external control. The latter is the subject of this paper, in which we study batch bioprocessing from a resource allocation perspective. Based on previous works, we propose a simple bacterial growth model accounting for the dynamics of the bioreactor and intracellular composition, and we analyze its asymptotic behavior and stability. Using optimization and optimal control theory, we study the production of biomass and metabolites of interest for infinite- and finite-time horizons. The resulting optimal control problems are studied using Pontryagin's Maximum Principle and numerical methods, and the solutions found are characterized by the presence of Fuller phenomenon (producing an infinite set of switching points occurring in a finite-time window) at the junctions with a second-order singular arc. The approach, inspired in biotechnological engineering, aims to shed light upon the role of cellular composition and resource allocation during batch processing and, at the same time, poses very interesting and challenging mathematical problems.

**Key words.** mathematical systems theory; nonlinear systems; mathematical cell model dynamics and control; industrial biotechnology; optimal control; bacterial resource allocation

MSC codes. 37N25, 49K15, 92C42

**1. Introduction.** The study of living microorganisms through resource allocation models has become increasingly relevant for its capacity to elucidate natural behaviors of microbia through very simple dynamical models [1, 3, 8, 12]. The core idea is to represent the distribution of cellular resources through optimal control strategies, based on the assumption that evolutionary processes have tuned these endogenous allocation strategies to attain nearly-optimal levels [14]. Numerous problems arise in this context, one of them being the optimal production of metabolites regulated by an external control capable of arresting bacterial growth [11]. Growth control has proven a key engineering method for several industrial applications, such as in food preservation, biofuel production, and in combating antibiotics resistance [10]. To this end, a resource allocation approach can help understand how to modify the naturally-evolved allocation strategies so as to efficiently produce such chemical compounds [6].

These biosynthetic strategies have been studied in different frameworks. The simplest case describes the interactions between intracellular proteins with minimal interplay with the environment [27, 5, 22]. The latter can be modelled by omitting the dynamics of the substrate in the medium, representing the case where bacterial exponential growth can be attained. Another relevant, more complex case is continuous bioreactors [23, 26], used extensively in industries and in cell biology research for its capacity to reach and maintain steady-state growth conditions. The latter is

<footnote>*Submitted to the editors 30th June 2022.</footnote>

<footnote>**Funding:** This work was partially supported by ANR project Maximic (ANR-17-CE40-0024-01), Ctrl-AB (ANR-20-CE45-0014), Inria IPL Cosy and Labex SIGNALIFE (ANR-11-LABX-0028-01).</footnote>

<footnote>†MISTEA, Université Montpellier, INRAE, Institut Agro, Montpellier, France (agustin.yabo@inrae.fr, www.agustinyabo.com.ar)</footnote>

<footnote>‡Université Côte d'Azur, CNRS, Inria, LJAD, France (jean-baptiste.caillau@univ-cotedazur.fr)</footnote>

<footnote>§Université Côte d'Azur, Inria, INRAE, CNRS, Sorbonne Université, Macbes Team, Sophia Antipolis, France (jean-luc.gouze@inria.fr)</footnote>
accomplished through an inflow of fresh medium rich in substrate and an outflow of the culture at the same volumetric flow rate, which produce a constant volume of the culture in the device. In that case, optimization studies are mostly oriented to reach such steady state in a cost-effective way. In fed-batch fermentation, the process starts with an initial volume of bacterial culture inside a bioreactor, which is progressively filled up through an inflow of rich medium, increasing the volume of the culture until it reaches a maximum level [24]. Once the maximum volume is attained, the culture evolves as a closed process, known in the field as batch processing. As no mass comes in or out of the device, the remainder of the nutrients in the medium are progressively consumed until the mass is entirely transformed into final products.

The latter is the subject of this paper, which tackles batch processing from a resource allocation perspective. The novelty of the approach lies in the nature of the model that—in addition to the physical and chemical laws found in classical bioreactor models—considers cellular composition, taking into account the intracellular components responsible for the main biological functions of bacteria. The problem has been first posed in [27], where a simpler mathematical model of resource allocation is studied through numerical optimal control. The study does not consider the dynamical aspects of the model, neither the theoretical specifics arising from the optimal control problem and its singular arcs. We extend these results from an analytical perspective—both for the dynamical analysis and the optimal control study—and including the case with no metabolite synthesis as a starting point, which has not been analyzed in previous works. Based on simpler bacterial growth models [7, 25] that do not consider the dynamics of the substrate in the medium, a coarse-grained self-replicator model is introduced, including a heterologous pathway for the production of a value-added chemical compound [27, 22]. Additionally, the main biological assumptions of the mechanistic bacterial model are revised, based on empirical studies of exponentially growing *E. coli* cultures [16]. Specifically, we consider a class of growth rate-independent proteins in the cellular composition that accounts for housekeeping proteins and non-active ribosomes, known to take up more than 50% of the cell [17]. The inclusion of this class of proteins in previous models has shown considerable improvement in the agreement between simulations and experimental data [25]. Using mass conservation laws related to the closeness of the bioprocess, it is possible to analyze the asymptotic behavior and stability of the dynamical system, showing that, for every possible allocation strategy, all component of the system are transformed either into proteins or into metabolites, a condition later defined as *Full depletion*. Then, two main studies are performed: the biomass maximization case, representing the natural objective of wild-type (i.e not modified) microbial cultures; and the metabolite maximization case, using the full bacterial model that includes the pathway for metabolite synthesis for industrial purposes. Both problems are analyzed in infinite time and in finite time, the latter stated as OCPs (Optimal Control Problems), which are investigated through the application of PMP (Pontryagin's Maximum Principle) [15]. While the finite-time case is suitable for representing bioprocesses with predetermined duration, the analysis of the infinite-time case becomes crucial in understanding the nature and asymptotic trend of the process. The solutions of the OCPs are characterized by the presence of Fuller's phenomenon [3], producing arcs composed of an infinite set of switching points (i.e. bangs) over a finite-time window. These optimal solutions follow a Fuller-singular-Fuller structure, similar to the one found in [25], described by a single second-order singular arc which is delimited by two Fuller's arcs at the beginning and at the end of the process. In particular, the solution of the biomass maximization case is thoroughly studied from an analytical
point of view, resulting in an explicit expression of the singular control in feedback
form. The results here presented are also confirmed by simulations obtained with
Bocop [18], an optimal control solver based on direct methods, and published in the
ct gallery<fnref n="1" /> in order to guarantee the reproducibility of the numerical results.

The paper is organized as follows: in Section 2, the dynamical model is pre-
sented, and its dynamical behavior is studied in Section 3. The biomass and product
maximization cases are introduced and investigated in Sections 4 and 5, respectively.
Finally, the results are discussed in Section 6.

## 2. Model definition.

**2.1. Self-replicator model.** We define a self-replicator model describing the dynamics of a microbial population growing inside a closed bioreactor. The bacterial culture has constant volume $V_e$, measured in liters. At the beginning of the experience, there is an initial mass of substrate S inside the bioreactor, that is gradually consumed by the bacterial population, and transformed into precursor metabolites P. These precursors are intermediate metabolites used to produce proteins—such as ribosomes and enzymes—responsible for specific cellular functions; and metabolites of interest X which are excreted from the cell. The proteins forming bacterial cells are divided into three classes M, R and Q, associated to the following cellular functions:

**Class M** Proteins of the metabolic machinery, responsible for the uptake of nutrients S from the medium, the production of precursor metabolites P, and the synthesis of metabolites of interest X.

**Class R** Proteins of the gene expression machinery (such as ribosomes) actively involved in protein biosynthesis (i.e. in the production of proteins of classes M, Q and R).

**Class Q** Growth rate-independent proteins, such as housekeeping proteins responsible for cell maintenance, and ribosomes not involved in protein synthesis [17].

From a biological perspective, the production of proteins M, R and Q is catalyzed by
ribosomal proteins R, and the absorption of S and synthesis of X are both catalyzed
by the metabolic proteins M. This catalytic effect is represented in Figure 1 through
dashed arrows. Intracellular proteins are produced at a synthesis rate $V_R$ measured
in grams per hour. The synthesis rates of proteins M, R and Q are $r_{max}(1-u)V_R$,
$r_{max}uV_R$ and $(1-r_{max})V_R$, respectively; where the parameter $r_{max}$ is a certain em-
pirical constant imposing a maximum threshold to the rate of production of proteins
M and R. The proportion of precursors dedicated to the production of growth rate-
independent proteins Q is fixed, while the balance between proteins M and R is
decided by the allocation control $u$. The latter is modelled through a time-varying
function $u(t) \in [0, 1]$, where $u = 0$ means no production of ribosomal proteins R, and
$u = 1$ means no production of metabolic proteins M. Depending on the objective to
be analyzed, the control $u$ can represent different mechanisms. First, it can account
for the natural allocation used by bacteria, as modelled in [7, 25], by assuming that
the native regulatory mechanisms of bacterial cells have been tuned by the natural
selection to maximize growth rate. On the other hand, it can represent the artifi-
cially modified allocation modelled in [27]. In a biotechnological setting, the latter
is accomplished by engineering a synthetic growth switch that allows to modify the
natural allocation through external compounds like IPTG<fnref n="2" />.

<footnote n="1">¹ct.gitlabpages.inria.fr/gallery/substrate/depletion.html</footnote>

<footnote n="2">²Isopropyl β-D-1-thiogalactopyranoside</footnote>
<figure><img src="image_2.png" /><figcaption>FIG. 1. Self-replicator model of bacterial growth representing the intracellular micro-chemical reactions behind nutrient uptake, cell growth and metabolite synthesis. Solid arrows represent flow of resources resulting from the microchemical reactions, while dashed arrows indicate a catalyzing effect (i.e. the presence of a protein accelerating the synthesis of another protein).</figcaption></figure>

2.2. Dynamical system. The dynamics of the self-replicator system are described by

$$
\begin{equation}
\begin{aligned}
\dot{S} &= -V_M, \\
\dot{P} &= V_M - V_X - V_R, \\
\dot{R} &= r_{\max} u V_R, \\
\dot{M} &= r_{\max} (1 - u) V_R, \\
\dot{Q} &= (1 - r_{\max}) V_R, \\
\dot{X} &= V_X,
\end{aligned}
\tag{2.1}
\end{equation}
$$

where the variables $S(t)$, $P(t)$, $R(t)$, $M(t)$, $Q(t)$ and $X(t)$ represent the masses (in 142 grams) of substrate, precursors metabolites, the gene expression machinery, the metabolic machinery, the growth rate-independent proteins and the metabolites of interest 143 at time $t$ measured in hours, respectively. $V_M(t)$, $V_R(t)$ and $V_X(t)$ are the reaction rates of the system (in grams per hour), and $u(t)$ is the allocation control previously 145 defined. We define the volume (in liters) of the bacterial population in the bacterial culture $V(t)$ as 146

$$ V \doteq \beta(M + R + Q), \tag{2.1} $$

where $\beta$ is a constant relating protein density and volume [2]. Definition (2.1) purposely neglects the mass of precursor metabolites $P(t)$, which greatly simplifies the computations. The latter assumption is based on the fact that most of the mass in 151 bacterial cells corresponds to proteins of classes M, R and Q, as confirmed in previous 152
155 studies [7]. This allows to define time-varying intracellular concentrations (in grams
156 per liter) with respect to this volume

$$
\begin{equation} \tag{2.2}
\begin{alignedat}{2}
& p &&\doteq \frac{P}{V}, \qquad r &&\doteq \frac{R}{V}, \qquad m &&\doteq \frac{M}{V}, \qquad q &&\doteq \frac{Q}{V}.
\end{alignedat}
\end{equation}
$$

159 Likewise, we define the extracellular concentrations related to the external volume

$$
\begin{equation}
\tag{2.3}
\begin{alignedat}{2}
s &= \frac{S}{V_e}, & x &= \frac{X}{V_e}.
\end{alignedat}
\end{equation}
$$

162 We define the relative synthesis rates involved in the processes as increasing functions
163 of the concentrations used in each reaction [17], and taking into account the catalytic
164 effect previously described

$$
\begin{equation}
v_M(s,m) \doteq \frac{V_M}{\mathcal{V}}, \quad v_R(p,r) \doteq \frac{V_R}{\mathcal{V}}, \quad v_X(p,m) \doteq \frac{V_X}{\mathcal{V}}.
\end{equation}
$$

167 From (2.1) and (2.2), we have that

$$
\begin{equation}
\tag{2.4}
m + r + q = \frac{1}{\beta},
\end{equation}
$$

170 which implies that the concentrations m, r and q cannot be bigger than 1/β. We
171 define the growth rate of the bacterial culture μ as

$$
\mu \doteq \frac{\dot{V}}{V} = \beta v_R(p, r).
$$

174 Then, the dynamical system can be expressed in terms of the concentrations as

$$
\left\{
\begin{aligned}
\dot{s} &= -v_M(s,m) \frac{V}{V_e}, \\
\dot{p} &= v_M(s,m) - v_X(p,m) - v_R(p,r)(\beta p + 1), \\
\dot{r} &= (r_{\max}u - \beta r)v_R(p,r), \\
\dot{m} &= (r_{\max}(1-u) - \beta m)v_R(p,r), \\
\dot{q} &= ((1-r_{\max}) - \beta q)v_R(p,r), \\
\dot{V} &= \beta v_R(p,r)V, \\
\dot{x} &= v_X(p,m) \frac{V}{V_e}.
\end{aligned}
\right.
$$

177 **2.3. Kinetics definition.** We model the kinetics of the system by supposing
178 that both the synthesis rates of precursors $v_M$ and metabolites $v_X$ are linear in the
179 concentration of metabolic proteins $m$, and the protein synthesis rate $v_R$ is linear in
180 the concentration of active ribosomal proteins $r$ [16]. Thus, they can be expressed as

$$
181 \qquad v_M(s,m) = w_M(s)m,
$$

$$
182 \qquad v_R(p, r) = w_R(p)r,
$$

$$
183 \qquad v_X(p, m) = \gamma w_R(p) m,
$$
where $\gamma > 0$ is a proportionality constant, which allows the metabolite synthesis rate to be expressed as $v_X(p, m) = \gamma v_R(p, r)m/r$. Such assumption implies that the bacterial cell has the same affinity to synthesize biomass and metabolites from the precursors, even if the reactions do not consume P in the same proportion. In the particular case of Michaelis-Menten kinetics, this feature is captured by the half-saturation constant [8]. The functions $w_I$ are assumed to have the following properties:

Hypothesis 2.1. Function $w_I(x) : \mathbb{R}_+ \to \mathbb{R}_+$ is

* Continuously differentiable w.r.t. $x$,

* Null at the origin: $w_I(0) = 0$,

* Strictly monotonically increasing: $w'_I(x) = \frac{\partial}{\partial x}w_I(x) > 0, \forall x \ge 0$,

* Strictly concave downwards: $w''_I(x) = \frac{\partial^2}{\partial x^2}w_I(x) < 0, \forall x \ge 0$,

* Upper bounded: $\lim_{x \to \infty} w_I(x) = k_I > 0$.

For numerical simulations, we resort to the particular case where the functions follow Michaelis-Menten kinetics. For that case, we define

$$w_R(p) \doteq k_R \frac{p}{K_R + p}, \quad w_X(p) \doteq k_X \frac{p}{K_X + p}, \quad w_M(s) \doteq k_M \frac{s}{K_M + s},$$

where the values of the constants $k_R, K_R, k_X, K_X, k_M$ and $K_M$ are based on the literature [7, 27]. For the general case introduced in Hypothesis 2.1, we define

$$k_R \doteq \lim_{p \to \infty} w_R(p), \quad k_X \doteq \lim_{p \to \infty} w_X(p), \quad k_M \doteq \lim_{s \to \infty} w_M(s).$$

**2.4. Mass fraction formulation and non-dimensionalization.** We define non-dimensional mass fractions

$$(2.5) \qquad \hat{s} \doteq \beta s, \quad \hat{p} \doteq \beta p, \quad \hat{r} \doteq \frac{\beta}{r_{\max}} r, \quad \hat{m} \doteq \frac{\beta}{r_{\max}} m, \quad \hat{q} \doteq \beta q, \quad \hat{x} \doteq \beta x,$$

where $\hat{r}$ and $\hat{m}$ are the mass fractions of the maximal ribosomal fraction $r_{\max}$. Then, given that the transcription of housekeeping proteins in bacterial cells is internally auto-regulated [20], and that the mass fraction of non-translating ribosomal proteins is constant [16], we assume that the mass fraction of growth rate-independent proteins $\hat{q}$ varies mildly compared to the remaining states, and thus we fix

$$(2.6) \qquad \hat{q} = 1 - r_{\max},$$

which, replacing in (2.4), yields

$$\hat{m} + \hat{r} = 1.$$

The latter implies that the metabolic fraction can be expressed in terms of the ribosomal fraction as $\hat{m} = 1 - \hat{r}$, and so the dynamical equation of $\hat{m}$ can be removed from the system. Additionally, we see that the quantity $\beta r$ represents the mass fraction of translating ribosomal proteins in the cell which, using (2.4) and (2.6), has bounds $[0, r_{\max}]$. Thus, its upper bound is given by the difference between the maximal total ribosomal mass fraction and the constant non-translating ribosomal mass fraction. In the literature [25], such values are empirically fixed to 0.5 and 0.07, respectively, and so the parameter $r_{\max}$ is here set to 0.43 for the numerical calculations. The biomass fraction of the bacterial culture is defined as

$$(2.7) \qquad \hat{V} \doteq \frac{V}{V_e}.$$
We define the non-dimensional time $\hat{t} \doteq k_R r_{\max} t$ and the non-dimensional functions

$$ \hat{w}_R(\hat{p}) = \frac{w_R(p)}{k_R}, \quad \hat{w}_X(\hat{p}) = \frac{w_X(p)}{k_R}, \quad \hat{w}_M(\hat{s}) = \frac{w_M(s)}{k_R} $$

so that $\lim_{\hat{p} \to \infty} \hat{w}_R(\hat{p}) = 1$. For the sake of simplicity, let us drop all hats from the current notation. Thus, the system becomes

$$ \begin{cases} \dot{s} = -w_M(s)(1-r)\mathcal{V}, \\ \dot{p} = w_M(s)(1-r) - \gamma w_R(p)(1-r) - w_R(p)r(p+1), \\ \dot{r} = (u-r)w_R(p)r, \\ \dot{\mathcal{V}} = w_R(p)r\mathcal{V}, \\ \dot{x} = \gamma w_R(p)(1-r)\mathcal{V}. \end{cases} \quad (S) $$

In this formulation, and using (2.2), (2.3), (2.4), (2.5) and (2.7), the total mass in the bioreactor can be expressed in terms of the concentrations as

$$ (2.8) \qquad S + P + M + R + Q + X = \frac{\mathcal{V}_e}{\beta}(s + (p+1)\mathcal{V} + x). $$

### 3. Model analysis.

**LEMMA 3.1.** *The set*

$$ \Gamma = \{(s, p, r, \mathcal{V}, x) \in \mathbb{R}^5 : s \ge 0, p \ge 0, 1 \ge r \ge 0, \mathcal{V} \ge 0, x \ge 0\} $$

is positively invariant for the initial value problem.

Proving Lemma 3.1 is standard and can be done by evaluating the vector field of (S) over the boundaries of $\Gamma$. Thus, we fix initial conditions

$$ \begin{cases} s(0) = s_0 > 0, \quad p(0) = p_0 > 0, \quad x(0) = 0, \\ r(0) = r_0 \in (0, 1), \quad \mathcal{V}(0) = \mathcal{V}_0 > 0, \end{cases} \quad (\text{IC}) $$

where the initial concentration of metabolites $x(0)$ is set to 0 to represent the fact that, at the beginning of the bioprocess, no metabolite has been produced. Some relations are immediate from the dynamics: as $\dot{s} \le 0$ and $\dot{\mathcal{V}} \ge 0$ for all $t$, we have

$$ (3.1) \qquad s(t) \le s_0, \quad \mathcal{V}(t) \ge \mathcal{V}_0, $$

representing the fact that the substrate can only be consumed (and not replenished), and the biomass can only grow.

**3.1. Total available mass.** As typically occurs in batch processes, there is neither inflow nor outflow of mass in the bioreactor, which is reflected in the dynamics of the system though a mass conservation law. We define the constant

$$ \Sigma \doteq s_0 + (p_0 + 1)\mathcal{V}_0, $$

representing the initial mass concentration in the system. It can be seen that the total mass concentration

$$ z \doteq s + (p+1)\mathcal{V} + x $$
262 is constant for all $t$ (as $\dot{z} = 0$). This means that

$$ (3.2) \qquad s + (p+1)\mathcal{V} + x = \Sigma, \qquad 263 $$

265 for all $t$. Thus, relation (2.8) and (3.2) show that the total mass in the system is 266 constant and equal to $\mathcal{V}_e\Sigma/\beta$. Variables $\mathcal{V}$ and $x$ are maximal when the remaining 267 variables are equal to 0, and so they are upper bounded. In particular, both $\mathcal{V}(t)$ and 268 $x(t)$ are decreasing w.r.t. $s(t)$ and $p(t)$. As neither $s$ nor $p$ can be negative, we have 269 that

$$ (3.3) \qquad \mathcal{V}(t) + x(t) = \Sigma, \qquad 270 $$

272 when $s(t) = p(t) = 0$. This condition means that all the available substrate and pre- 273 cursor metabolites have been depleted and transformed into biomass and metabolites, 274 which is intuitively what one would expect from system (S) for $t$ sufficiently large. 275 Additionally, using (3.1) and (3.2), we can obtain the following result.

276 PROPOSITION 3.2. $\mathcal{V}(t) \in [\mathcal{V}_0, \Sigma]$, $x(t) \in [0, \Sigma - \mathcal{V}_0]$ and $p(t) \in [0, p^+]$ for all $t$, 277 with $p^+ = \Sigma/\mathcal{V}_0 - 1$.

278 3.2. Infinite-time full depletion. Dynamics (S) shows that, under initial con- 279 ditions (IC), $s(t)$ and $p(t)$ can only vanish asymptotically, that is, when $t \to \infty$. The 280 latter can be proved by seeing that the derivatives of $s$ and $p$ can be bounded by

$$ \dot{s} \ge -w_M(s)\Sigma, \quad \dot{p} \ge -w_R(p)p^+(p^+ + 1 + \gamma), \qquad 281 $$

283 which means that, at worst, $s$ and $p$ decay exponentially (as functions $w_i(x)$ can be 284 upper bounded by linear functions $w_i(x) \le c_i x$), so that $s(t) = p(t) = 0$ cannot be 285 reached in finite time. Accordingly, we define the depletion of $s$ and $p$ in an infinite- 286 time horizon.

287 DEFINITION 3.3. System (S) achieves Full depletion when all the substrate and 288 the precursors are asymptotically depleted, i.e.

$$ (\text{Full depletion}) \qquad \lim_{t \to \infty} s(t) = \lim_{t \to \infty} p(t) = 0, \qquad 289 $$

291 3.3. Asymptotic behaviour. Now, we study the system dynamics for an infi- 292 nite time $t \to \infty$. First, the case with a constant allocation $u(t) = u^* \in (0, 1)$ for all 293 $t$ is analyzed, and then an extension to a general allocation function is proposed.

294 3.3.1. Constant allocation $u^*$.

295 THEOREM 3.4. For any trajectory of system (S) with initial conditions (IC) and 296 constant allocation $u(t) = u^*$, it follows that

$$ (3.4) \qquad (u^* - r(t))\mathcal{V}(t) = (u^* - r_0)\mathcal{V}_0. \qquad 297 $$

Proof. Under a constant allocation $u(t) = u^*$, the dynamics of $r$ becomes

$$ \dot{r} = (u^* - r)w_R(p)r. $$

Using dynamics (S), it is possible to see that both the total mass of proteins $R = r\mathcal{V}$ and the quantity $U = u^*\mathcal{V}$ have the same derivative

$$ \dot{R} = \dot{U} = u^*w_R(p)r\mathcal{V}, $$

299 which means that the difference of these two $R_u = U - R$ should be constant (as 300 $\dot{R}_u = 0$), which yields (3.4). □
### 3.3.2. General allocation $u(t)$. Due to the boundedness of $V$ stated in Lemma 3.2, and the relation between $V$ and $r$ shown in (3.2), we can see that any constant control $u^*$ yields a bounded ribosomal fraction $r$. We extend this notion to any function $u(t)$.

**LEMMA 3.5.** For any trajectory of system (S) with initial conditions (IC) and any control $u(t)$, the ribosomal concentration has bounds $r(t) \in [r^-, r^+]$ for all $t$, with

$$r^- \doteq r_0 \frac{V_0}{\Sigma} > 0, \quad r^+ \doteq 1 - (1 - r_0) \frac{V_0}{\Sigma} < 1.$$

*Proof.* Let us extend system (S) by defining variables $r_{\text{low}}(t)$ and $r_{\text{up}}(t)$ with dynamics

$$\begin{aligned} \dot{r}_{\text{low}} &= -r_{\text{low}}w_R(p)r \le 0, & \dot{r}_{\text{up}} &= (1 - r_{\text{up}})w_R(p)r \ge 0, \\ r_{\text{low}}(0) &= r_0, & r_{\text{up}}(0) &= r_0, \end{aligned}$$

which correspond to the dynamics of $r$ with $u=0$ and $u=1$ respectively, and which satisfy

$$r_{\text{low}}(t) \le r(t) \le r_{\text{up}}(t)$$

for all $t$. The latter can be easily proved by showing that the time-varying differences

$$\Delta_{\text{low}}(t) = r(t) - r_{\text{low}}(t), \quad \Delta_{\text{up}}(t) = r_{\text{up}}(t) - r(t)$$

with dynamics

$$\dot{\Delta}_{\text{low}} = (u - \Delta_{\text{low}})w_R(p)r, \quad \dot{\Delta}_{\text{up}} = (1 - u - \Delta_{\text{up}})w_R(p)r$$

are always non-negative: they satisfy $\Delta_{\text{low}}(0) = \Delta_{\text{up}}(0) = 0$ and are repulsive or (at worst) invariant at 0. Then, based on the same principle used to obtain (3.4), we define the quantities $R_{\text{low}} = r_{\text{low}}V$ and $R_{\text{up}} = (1 - r_{\text{up}})V$ which are constant (as $\dot{R}_{\text{low}} = \dot{R}_{\text{up}} = 0$), and so

$$r_{\text{low}}(t) = r_0 \frac{V_0}{V(t)}, \quad r_{\text{up}}(t) = 1 - (1 - r_0) \frac{V_0}{V(t)},$$

for all $t$. As $V_0 \le V(t) \le \Sigma$ for all $t$, we have

$$r_{\text{low}}(t) \in \left[ r_0 \frac{V_0}{\Sigma}, r_0 \right], \quad r_{\text{up}}(t) \in \left[ r_0, 1 - (1 - r_0) \frac{V_0}{\Sigma} \right]$$

which shows that $r^- \le r(t) \le r^+$ for all $t$. □

**LEMMA 3.5.** States that, for any control $u(t)$, the ribosomal concentration never reaches the bounds $r=0$ and $r=1$, and thus neither the substrate intake nor the protein synthesis is arrested. Using this fact, it can be proved that any control $u(t)$ produces (Full depletion).

**THEOREM 3.6.** Any trajectory of system (S) with initial conditions (IC) and any control $u(t)$ achieves (Full depletion) when $t \to \infty$.

*Proof.* Using Lemma 3.5, it is easy to see that

$$\dot{s} \le -w_M(s)(1 - r^+)V_0,$$

which means that $s(t)$ converges to 0 as $t \to \infty$. Then, this means that

$$\dot{p} \le -\gamma w_R(p)(1 - r^+) - w_R(p)r^-,$$

and so $p(t)$ also converges to 0 as $t \to \infty$. □
**4. The biomass maximization case.** In this section, we write the problem of maximizing the biomass both for infinite time and finite time in terms of the allocation parameter $u$. The latter is a mathematical representation of the naturally-evolved resource allocation strategy used by bacteria in nature. Indeed, in biology it is very often assumed that bacteria during exponential growth allocate their internal resources to maximize their growth rate, thus maximizing long-term biomass production [7]. For this particular problem, we assume that no metabolite is produced, as the pathway responsible for its production is artificially engineered, and thus not present in wild-type bacteria. This is simply modeled through $\gamma = 0$. The resulting Wild-Type Bacterial Model is

$$ (\text{WTB-M}) \quad \left\{ \begin{aligned} \dot{s} &= -w_M(s)(1-r)\mathcal{V}, \\ \dot{p} &= w_M(s)(1-r) - w_R(p)r(p+1), \\ \dot{r} &= (u-r)w_R(p)r, \\ \dot{\mathcal{V}} &= w_R(p)r\mathcal{V}, \end{aligned} \right. $$

**4.1. Infinite-time problem.**

**4.1.1. Problem formulation.** We first write the biomass maximization problem for an infinite-time horizon, a non-realistic scenario that can provide valuable insight into the finite-time process. Indeed, in this section, we show that the maximum attainable performance can only be achieved in infinite-time processes. The problem can be expressed as

$$ \max_{u(t)} \lim_{t \to \infty} \mathcal{V}(t). $$

Since $\mathcal{V} \in [\mathcal{V}_0, \Sigma]$, applying (Full depletion) in (3.2) yields the condition

$$ \lim_{t \to \infty} \mathcal{V}(t) = \Sigma. $$

meaning that, in infinite time, the biomass is maximized for every control $u(t)$. We formalize the latter in the following theorem.

**THEOREM 4.1.** For any trajectory of system (WTB-M) with initial conditions (IC) and any control $u(t)$, the volume $\mathcal{V}(t) \to \max \mathcal{V}(t) = \Sigma$ as $t \to \infty$.

As a consequence, using Theorem 3.4, we have the following result for constant allocations.

**COROLLARY 4.2.** For any trajectory of system (WTB-M) with initial conditions (IC) and constant control $u(t) = u^*$,

$$ \lim_{t \to \infty} r(t) = u^* - (u^* - r_0) \frac{\mathcal{V}_0}{\Sigma} $$

These results are illustrated by the numerical simulations shown in next section.

**4.1.2. Numerical simulations.** Examples of trajectories confirming the analytical results are shown in Figure 2 and Figure 3, where we see that the system approaches (Full depletion) asymptotically in every case, thus approaching the maximal biomass value $\mathcal{V}(t) = \Sigma$. Figure 2 shows the resulting trajectories associated to
373 the same initial conditions, when varying the allocation parameter $u$. On the other hand, Figure 3 illustrates the trajectories for different values of $r_0$. Indeed, as $\Sigma$ does 374 not depend on the resource allocation strategy, all the available mass is transformed 375 into biomass independently of the values of $r_0$ and $u(t)$. 376

<figure><img src="image_3.png" /><figcaption>FIG. 2. Simulation of (WTB-M) with initial conditions $s_0 = 0.3$, $p_0 = 0.001$, $r_0 = 0.8$, $V_0 = 0.003$, fixed final time $t_f = 50$ and different allocation functions $u$.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>FIG. 3. Simulation of (WTB-M) with initial conditions $s_0 = 0.3$, $p_0 = 0.001$, $V_0 = 0.003$, $u = 0.5$, fixed final time $t_f = 50$ and different values of $r_0$.</figcaption></figure>

## 4.2. Finite-time problem.
378

379

380

**4.2.1. Problem formulation.** For the Biomass Maximization problem at final time $t_f$, we write the OCP maximizing the final bacterial volume $V(t_f)$ with initial conditions (IC):

381

$$ \left\{ \begin{array}{ll} \text{maximize} & V(t_f), \\ \text{subject to} & \text{dynamics of (WTB-M),} \\ & \text{initial conditions (IC),} \\ & u(\cdot) \in \mathcal{U}. \end{array} \right. $$

382

383

384

385

386

387

388

389

390

For this class of optimal control problem, with no terminal constraints, there are no controllability issues. Additionally, the dynamics is affine in the control, with the latter included in a compact and convex set (a closed interval), and it can be checked that every finite-time trajectory remains bounded. Thus, existence of a solution is guaranteed by Filippov's theorem [1]. Then, for a problem (BM-OCP) with state $\varphi \in \mathbb{R}^n$, PMP ensures that there exists an absolutely continuous mapping $\lambda(\cdot) : [0, t_f] \to \mathbb{R}^n$ such that the extremal $(\varphi, \lambda, u)$ satisfies the generalized Hamiltonian system

391

$$ \left\{ \begin{array}{ll} \dot{\varphi} = \frac{\partial}{\partial \lambda} H(\varphi, \lambda, u), \\ \dot{\lambda} = -\frac{\partial}{\partial \varphi} H(\varphi, \lambda, u), \\ H(\varphi, \lambda, u) = \max_{u \in [0,1]} H(\varphi, \lambda, u), \end{array} \right. $$

392

393

394

for almost every $t \in [0, t_f]$. We define the adjoint states for this particular case as $\lambda = (\lambda_s, \lambda_p, \lambda_r, \lambda_\nu)$, and we write the Hamiltonian

395

396

$$ H = -w_M(s)(1-r)V\lambda_s + (w_M(s)(1-r) - w_R(p)r(p+1))\lambda_p + w_R(p)rV\lambda_\nu \\ + (u-r)w_R(p)r\lambda_r, $$

398

and the adjoint system as

399

$$ \left\{ \begin{array}{ll} \dot{\lambda}_s = w'_M(s)(1-r) (V\lambda_s - \lambda_p), \\ \dot{\lambda}_p = (w'_R(p)r(p+1) + w_R(p)r)\lambda_p - w'_R(p)rV\lambda_\nu - (u-r)w'_R(p)r\lambda_r, \\ \dot{\lambda}_r = -w_M(s)(V\lambda_s - \lambda_p) + w_R(p)((p+1)\lambda_p - V\lambda_\nu) - (u-2r)w_R(p)\lambda_r, \\ \dot{\lambda}_\nu = w_M(s)(1-r)\lambda_s - w_R(p)r\lambda_\nu. \end{array} \right. $$

400

401

402

403

404

405

Given that there are no terminal conditions on the state, the transversality conditions for the adjoint state are $\lambda(t_f) = (0, 0, 0, 1)$. Note that $\lambda_\nu(t_f) = 1$ comes from the fact that the cost function is $V(t_f)$ and there are no terminal conditions on the other states (this prevents the so-called abnormal extremals with $\lambda_\nu(t_f) = 0$). Since the (WTB-M) dynamics is single-input and control-affine,

406

$$ \dot{\varphi} = F_0(\varphi) + uF_1(\varphi) $$
with obvious definitions for the vector fields $F_0, F_1$, the Hamiltonian writes $H = H_0 + uH_1$. The Hamiltonian lifts $H_i := \langle \lambda, F_i \rangle$, $i = 0, 1$, of the two vector fields are

$$ 
\begin{aligned}
H_0 &= -w_M(s)(1-r)\mathcal{V}\lambda_s + (w_M(s)(1-r) - w_R(p)r(p+1))\lambda_p \\
    &\quad + w_R(p)r\mathcal{V}\lambda_\nu - w_R(p)r^2\lambda_r, \\
H_1 &= w_R(p)r\lambda_r.
\end{aligned}
$$

The constrained optimal control $u$ should maximize the Hamiltonian, so the solution of (BM-OCP) is

$$ (4.1) \qquad u(t) = \begin{cases} 0 & \text{if } H_1 < 0, \\ 1 & \text{if } H_1 > 0, \end{cases} $$

while $u(t) = u_s(t)$ is called singular whenever $H_1$ vanishes on a whole subinterval of $[0, t_f]$. This tells that an optimal control is a (possibly very complicated) concatenation of bangs ($u = 0$ and $u = 1$) and singular arcs, depending on the sign of the switching function $H_1$. We see that, at final time, the dynamics of $\lambda_r$ becomes

$$ \dot{\lambda}_r(t_f) = w_R(p(t_f))\mathcal{V}(t_f)\lambda_0 < 0, $$

which, using the fact that $\lambda_r(t_f) = 0$, implies that $\lambda_r > 0$ for a period $[t_f - \varepsilon, t_f]$, and thus $H_1 > 0$ for a period $[t_f - \varepsilon, t_f]$. As the control $u$ should maximize the Hamiltonian, we have the following result.

**LEMMA 4.3.** *There exists $\varepsilon$ such that the final bang arc of the optimal control solution of (BM-OCP) corresponds to $u(t) = 0$ for the time interval $[t_f - \varepsilon, t_f]$.*

A more detailed analysis of the optimal control solution can be done by studying the behavior of singular extremals. The latter is key in describing the structure of the optimal control, as it is typically associated to intermediate values of the control $u$ (i.e. non-bang arcs). In the general case application of PMP, the singular control can be expressed as a function of the state and the adjoint state $u_s(t) = f(\varphi, \lambda)$, where the explicit expression of $f$ can be obtained by successively differentiating the switching function $H_1$ until the singular control can be solved for. In the next section, we show for our problem that the singular optimal control is of order two, and that it can be expressed in feedback form as $u_s(t) = u(\varphi(t))$. In control systems design, the latter is a particular case that allows for a straightforward closed-loop implementation of the optimal control law, and that can provide further insight on the nature of the optimal trajectories.

**4.2.2. Singular arcs.** A singular arc occurs when $H_1$ vanishes (as well as its successive derivatives w.r.t. time) on a subinterval $[t_1, t_2] \subset [0, t_f]$, and so

$$ H_1 = w_R(p)r\lambda_r = 0, $$

As $r$ is bounded and $p$ cannot vanish in finite time, the condition becomes

$$ (4.2) \qquad \lambda_r = 0. $$

We differentiate $H_1$, we evaluate the expression in (4.2), and we get

$$ \dot{H}_1 = w_R(p)r \left( -w_M(s)(\mathcal{V}\lambda_s - \lambda_p) + w_R(p)((p+1)\lambda_p - \mathcal{V}\lambda_\nu) \right) = 0. $$
450 Then, the Hamiltonian can be expressed as

$$451 \quad (4.3) \qquad H = -w_M(s)(\mathcal{V}\lambda_s - \lambda_p) - r\lambda_r + (u-r)H_1 = c$$

on the interval $[t_1, t_2]$, where $c$ is a positive constant that, due to the constancy of the Hamiltonian, is equal to

$$c \doteq H(t_f) = w_R(p(t_f))r(t_f)\mathcal{V}(t_f) > 0.$$

453 Then, (4.3) implies that, over a singular arc,

$$454 \quad (4.4) \qquad w_M(s)(\mathcal{V}\lambda_s - \lambda_p) = w_R(p)((p+1)\lambda_p - \mathcal{V}\lambda_\nu) = -c.$$

456 Differentiating $\dot{H}_1$, and evaluating over $H_1 = \dot{H}_1 = 0$, yields

$$457 \quad \ddot{H}_1 = w_R(p)r(w'_R(p)((p+1)\lambda_p - \mathcal{V}\lambda_\nu)(w_M(s)(1-r) - w_R(p)r(p+1)) \\ 458 \qquad + w_R(p)(w_M(s)(1-r) - w_R(p)r(p+1))\lambda_p \\ 459 \qquad + w_R(p)(p+1)((w'_R(p)r(p+1) + w_R(p)r)\lambda_p - w'_R(p)r\mathcal{V}\lambda_\nu) \\ 460 \qquad - w_R^2(p)r\mathcal{V}\lambda_\nu - w_R(p)\mathcal{V}(w_M(s)(1-r)\lambda_s - w_R(p)r\lambda_\nu)) = 0.$$

462 By simplifying terms and replacing with (4.4), we obtain

$$463 \quad \ddot{H}_1 = cr(1-r)(w_R^2(p) - w'_R(p)w_M(s)) = 0.$$

465 which implies that

$$466 \quad (4.5) \qquad \frac{w_R^2(p)}{w_M(s)w'_R(p)} = 1.$$

468 The fact that $u$ does not appear in $\ddot{H}_1$ shows that any singular arc is at least of 469 local order two, so that additional derivatives should be calculated in order to retrieve 470 an explicit expression of the optimal control. Here, some precisions are in order, 471 and we may first recall that the computation can also be performed in terms of 472 Poisson brackets<fnref n="3" /> since derivating along an extremal amounts to bracketing with the 473 Hamiltonian. In particular,

$$474 \quad \dot{H}_1 = \{H_0 + uH_1, H_1\}, \\ 475 \quad = \{H_0, H_1\} =: H_{01}.$$

477 Iterating, and with obvious notations ($H_{001} := \{H_0, \{H_0, H_1\}\}$, etc.), one obtains

$$478 \quad 0 = \dot{H}_{01} = H_{001} + uH_{101}.$$

479 (Let us recall that $H_{01}$ is also equal to the Hamiltonian lift of the Lie bracket 480 $[F_0, F_1] =: F_{01}$, that $H_{001}$ is the lift of $[F_0, [F_0, F_1]] =: F_{001}$, etc.) The previous 481 computation shows that $H_{101}$ is zero on the subset $\{H_1 = H_{01} = 0\}$ of the cotangent 482 space. These two relations have indeed been used during the computation, while an

<footnote n="3">³In coordinates, if $fg$ and $g$ are two scalar valued functions of $(\varphi, \lambda) \in \mathbb{R}^{2n}$, $\{f, g\} = \sum_{i=1}^{n}(\partial f/\partial \lambda_i)(\partial g/\partial \varphi_i) - (\partial f/\partial \varphi_i)(\partial g/\partial \lambda_i)$.</footnote>
explicit evaluation<fnref n="4" /> allows to verify that the bracket $H_{101}$ is not identically zero on the whole cotangent space. In such a situation, the local (not intrinsic) order is said to be at least two; that is at least two more differentiations wrt. time are required to retrieve the singular control. We are actually in the following case (see also [4] for a similar analysis):

**PROPOSITION 4.4.** If the Lie bracket $F_{101}$ belongs to the span of $F_1$ and $F_{01}$, then singular extremals must be of (local) order at least two.

*Proof.* By assumption, if for some $\varphi$ a covector $\lambda$ is orthogonal to $F_1$ and $F_{01}$ at $\varphi$, it is also orthogonal to $F_{101}$ at this point. So, along a singular extremal that must belong to $\{H_1 = H_{01} = 0\}$, one has

$$0 \equiv H_{001} + u_s H_{101}$$

with $H_{101}$ also vanishing. As a result, $0 \equiv H_{001}$ along the singular, and one can differentiate again:

$$0 \equiv H_{0001} + u_s H_{1001}.$$

Now, by Leibniz rule

$$H_{1001} = \{H_1, \{H_0, H_{01}\}\} = \underbrace{\{-H_{01}, H_{01}\}}_{=0} + \{H_0, H_{101}\}$$

and there exist smooth functions $a$ and $b$ of $\varphi$ such that $F_{101} = aF_1 + bF_{01}$ (and similarly for the associated Hamiltonian lifts). By Leibniz rule again, $H_{1001}$ must vanish when $H_1 = H_{01} = H_{001} = 0$ as

$$\{H_0, aH_1 + bH_{01}\} = \{H_0, a\}H_1 + aH_{01} + \{H_0, b\}H_{01} + bH_{001},$$

so $H_{0001} \equiv 0$. So one has to differentiate at least once more to retrieve the control. □

It should moreover be noted that $\ddot{H}_1$ depends only on the state—and not on the adjoint state—which implies that its successive derivatives also depend only on the state, as the adjoint state does not appear in system (WTB-M). Additionally, and based on Hypothesis (2.1), the function $w_M(s)$ is invertible, which means that $s$ can be expressed in terms of $p$ through equation (4.5). Once again, we differentiate $\ddot{H}_1$, we evaluate over $H_1 = \dot{H}_1 = \ddot{H}_1 = 0$, and we get

$$\begin{aligned} \dot{\ddot{H}}_1 &= c(1-r)r(2w_R(p)w'_R(p)(w_M(s)(1-r) - w_R(p)r(p+1)) \\ &\quad - w''_R(p)w_M(s)(w_M(s)(1-r) - w_R(p)r(p+1)) \\ &\quad + w'_R(p)w'_M(s)w_M(s)(1-r)\mathcal{V}) = 0. \end{aligned}$$

By rearranging the expression, we can express

$$\dot{\ddot{H}}_1 = c(1-r)r(\omega_0(p, \mathcal{V}) - \omega_1(p, \mathcal{V})r) = 0,$$

<footnote n="4">⁴Take for instance $k_R = 1.1$, $k_M = 1.2$, $K_R = 1.3$, $K_M = 1.4$, $\varphi = \lambda = (1, 1, 1, 1)$ and check that $H_{101}(\varphi, \lambda) \neq 0$.</footnote>
517 with

$$518 \qquad \omega_0(p, \mathcal{V}) = w_M(s) (2w_R(p)w'_R(p) - w''_R(p)w_M(s) + w'_R(p)w'_M(s)\mathcal{V}) > 0,$$

$$519 \qquad \omega_1(p, \mathcal{V}) = \omega_0(p, \mathcal{V}) + w_R(p)(p+1)(2w'_R(p)w_R(p) - w''_R(p)w_M(s)) > 0,$$

521 where the positivity of functions $\omega_i$ is a consequence of Hypothesis 2.1. The latter shows that, along the singular arc, $r$ can be expressed in terms of $p$ and $\mathcal{V}$ (as, using (4.5), $s$ can be expressed in terms of $p$) as

$$524 \qquad r = \frac{\omega_0(p, \mathcal{V})}{\omega_1(p, \mathcal{V})}.$$

526 Then, computing the next derivative and evaluating over the obtained conditions 527 yields

$$528 \qquad (4.6) \qquad \ddot{H}_1 = c(1-r)r(\dot{\omega}_0(p, \mathcal{V}) - \dot{\omega}_1(p, \mathcal{V})r - \omega_1(p, \mathcal{V})(u-r)w_R(p)r) = 0.$$

530 We can see that the factor of $u$ in expression (4.6) satisfies

$$531 \qquad (4.7) \qquad \frac{\partial}{\partial u} \ddot{H}_1 = -c(1-r)w_R(p)r^2 \omega_1(p, s, \mathcal{V}) < 0$$

533 as $\omega_1 > 0$ for all $t$, which leads to the following result.

534 THEOREM 4.5. *The singular optimal control is exactly of order two, and it can be 535 expressed in feedback form, $u = u(p, \mathcal{V})$.*

536 *Proof.* Since it has been already proven that the singular arc is at least of order 537 two, it suffices to prove that it cannot be of greater order. This can be done by 538 observing that (4.7) cannot vanish. Thus, the singular arc is exactly of order two. 539 Additionally, solving for $u$ in the same expression yields

$$540 \qquad u_s(p, \mathcal{V}) = \frac{\dot{\omega}_0(p, \mathcal{V}) - \dot{\omega}_1(p, \mathcal{V})r}{\omega_1(p, \mathcal{V})w_R(p)r} + r.$$

542 which shows that the control $u$ can be expressed as a function of the state variables.□

543 We note that the generalized Legendre-Clebsh condition, necessary for optimality of 544 the singular arc, is fulfilled in strict form by virtue of (4.7):

$$545 \qquad (-1)^k \frac{\partial}{\partial u} \left( \frac{d^{2k}}{dt^{2k}} H_1 \right) = \frac{\partial}{\partial u} \ddot{H}_1 < 0.$$

547 **4.2.3. Numerical simulations.** The optimal trajectories were computed with 548 Bocop [18], which solves the OCP through a direct method. The time discretization 549 algorithm used is Lobato IIIC (implicit, 4-stage, order 6) with 2000 time steps. Figures 550 4 and 5 show optimal trajectories for the same set of initial conditions and different 551 values of $t_f$. Using the mass conservation law (3.2), the quantities are represented in 552 the plots as fractions of the total mass in the bioreactor $V_e\Sigma/\beta$. The optimal control 553 $u$ is characterized by the presence of chattering after and before the singular arc, as 554 expected in singular arcs of order two. From a biological point of view, both allocation 555 strategies prioritize the synthesis of proteins of the metabolic machinery M (red in 556 both Figures): the singular arc takes rather small values, and a large proportion of 557 the optimal control corresponds to a bang arc $u = 0$ at the end of the bioprocess.
The latter strategy promotes nutrient uptake, which results in a faster depletion of the substrate. It is interesting to note that, in opposition to previous results in the literature [7, 25], the presence of the turnpike properties [19] is not assured. Intuitively, the turnpike phenomenon would cause the time interval corresponding to the singular arc to increase as the final time $t_f$ increases. However, a quick comparison between Figures 4 and 5 shows that this is not the case, as the increase of $t_f$ only produces a larger 0-bang arc at the end of the bioprocess. This suggests that the duration of the initial phase dedicating a fraction of the resources to ribosomal proteins is fixed and independent of the duration of the batch process. Figure 6 also shows an optimal trajectory, but with a different initial ribosomal concentration. According to multiple simulations, this change shows an impact on the initial Fuller arc, that becomes perceptively larger, but has no major effects on the remaining of the process. Additionally, we note that the concentration of precursor metabolites in the bioreactor remains negligible in comparison with the other quantities, a result that matches the biological assumptions done in the modelling section.

<figure><img src="image_5.png" /><figcaption>FIG. 4. Numerical simulation of (BM-OCP) with initial conditions $s_0 = 0.1$, $p_0 = 0.001$, $r_0 = 0.1$, $V_0 = 0.003$ and $t_f = 30$. Quantities in the right plot are shown as fractions of the total mass in the bioreactor $V_e\Sigma/\beta$. The final volume $V(t_f)$ is at 95% of $V_e\Sigma/\beta$.</figcaption></figure>

Figure 7 illustrates optimal trajectories in the sp-plane for different final times (20, 25, 30 and 40). Each trajectory approaches the singular curve $\dot{H}_1 = 0$ given by expression (4.5) (obtained from the singular surface) through a Fuller arc, slides along it during a certain time interval, and then follows a trajectory obtained from the $u = 0$ arc that approaches asymptotically (Full depletion). Naturally, the longer the simulation of the process, the closer the final state to (Full depletion). These results also confirm the observations previously done: the duration in time of the singular arc is not directly related to the final time $t_f$. In fact, all the singular arcs start at approximately the same time instant and finish around $t = 18$. The independence of the initial phase from the duration of the bioprocess is coherent with the fact that bacteria allocate resources in terms of their cellular composition and the environment [7] (which, in this case, is described by the concentration of substrate in the medium), independently of the man-made notion of duration of the bioprocess. It is also noteworthy that, while the processes exit the singular surface at similar times, the trajectories differ significantly as they exit the singular arc from different initial conditions (not only in the $(s,p)$ plane but also in the original $\mathbb{R}^4$ space).
<figure><img src="image_6.png" /><figcaption>FIG. 5. Numerical simulation of (BM-OCP) with initial conditions $s_0 = 0.1$, $p_0 = 0.001$, $r_0 = 0.1$, $V_0 = 0.003$ and $t_f = 40$. Quantities in the right plot are shown as fractions of the total mass in the bioreactor $V_e\Sigma/\beta$. The final volume $V(t_f)$ is at 99.8% of $V_e\Sigma/\beta$.</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>FIG. 6. Numerical simulation of (BM-OCP) with initial conditions $s_0 = 0.1$, $p_0 = 0.001$, $r_0 = 0.3$, $V_0 = 0.003$ and $t_f = 40$. Quantities in the right plot are shown as fractions of the total mass in the bioreactor $V_e\Sigma/\beta$. The final volume $V(t_f)$ is at 99.8% of $V_e\Sigma/\beta$.</figcaption></figure>

### 4.2.4. Alternative approach: prescribed performance in minimum time.

As confirmed by the finite-time case studied in the last section, the associated optimal control problem with fixed final time yields a final volume $V(t_f)$ that can be viewed as a fraction (between 0 and 1) of the total mass concentration in the system $\Sigma$. Indeed, Theorem 4.1 showed that $V$ can reach its maximum value $\Sigma$ only when $t$ goes to infinity. Thus, (BM-OCP) can be reformulated to achieve a minimal-time transfer between an initial state (IC) (with biomass $V(0) = V_0$) and a final state with terminal constraint

589

590

591

592

593

594

595

596

597

(TC)

$$V(t_f) = \eta\Sigma,$$
<figure><img src="image_8.png" /><figcaption>FIG. 7. Numerical simulation of (BM-OCP) showing different trajectories in the sp-plane. Initial conditions are set to $s_0 = 0.1$, $p_0 = 0.001$, $r_0 = 0.1$, $V_0 = 0.003$. In all cases, the state approaches the singular curve $\ddot{H}_1 = 0$ and slides along it.</figcaption></figure>

for a certain performance parameter $\eta \in [\eta_{\min}, 1)$, where $\eta_{\min} \doteq V_0/\Sigma$. The reformu- lated minimal-time OCP with Prescribed Performance writes

$$ \text{(PP-OCP)} \quad \left\{ \begin{array}{ll} \text{minimize} & t_f, \\ \text{subject to} & \text{dynamics of (WTB-M)}, \\ & \text{initial conditions (IC)}, \\ & \text{terminal constraints (TC)}, \\ & u(\cdot) \in \mathcal{U}. \end{array} \right. $$

A natural question arising from OCPs with terminal constraints is the existence of a solution. In this case, Theorem 4.1 guarantees that any final volume $V(t_f) = \eta\Sigma$ can be reached in finite time, as long as $\eta \in [\eta_{\min}, 1)$. The latter ensures the existence of the solution for any $\eta \in [\eta_{\min}, 1)$. The study of the solutions of (PP-OCP) can be performed through an analogous PMP approach, with the difference that the Hamiltonian is null for every $t \in [0, t_f]$ due to the free final time $t_f$, and that there is no terminal constraint on $\lambda_v$ (i.e. $\lambda_v(t_f)$ is free). Given the similarity with the previously analyzed case, the computations of such optimal control solution are not explicated here. A numerical solution of the problem is shown in Figure 8.

**5. The product maximization case.** As done in the previous section, we approach the product maximization objective in infinite time and finite time using the full model (S) where $\gamma \in \mathbb{R}^+$.

**5.1. Infinite-time problem.** The problem of maximizing the product concen- tration at infinite time is given by the expression

$$ \max \lim_{u^* \to \infty} x(t), $$
<figure><img src="image_9.png" /><figcaption>FIG. 8. Numerical simulation of (PP-OCP) with initial conditions $s_0 = 0.1$, $p_0 = 0.003$, $r_0 = 0.1$ and $V_0 = 0.003$. The performance parameter is fixed to $\eta = 0.8$, which is achieved in $t_f = 25.1$. Quantities in the right plot are shown as fractions of the total mass in the bioreactor $V_e\Sigma/\beta$.</figcaption></figure>

which, using (3.3), can be rewritten as

$$ \min \lim_{u^* \to \infty} V(t) $$

indicating that maximizing the metabolite concentration at infinite time equates to minimizing the biomass. While, in the previous section, the conditions (Full depletion) and (3.2) were sufficient to determine the asymptotic behavior of the system, the presence of $x$ in this particular problem does not allow a similar resolution. An alternative approach, in order to better understand the role of cellular composition in the final objective, is to study a simplified version of the problem assuming cellular composition is at steady state (a common property of bacterial cells during exponential growth). We then propose a reduced version of the dynamical system by fixing the ribosomal concentration to a constant value $r^*$, which reduces the dimension of the model by one. Thus, the metabolite maximization problem is solved in terms of the parameter $r^*$, which represents a simpler analysis that can potentially provide an insight into the original optimization problem.

**5.1.1. Constant ribosomal concentration.** The system with Constant Ribosomal Concentration $r^*$ writes

$$ \begin{cases} \dot{s} = -w_M(s)(1-r^*)V, \\ \dot{p} = w_M(s)(1-r^*) - \gamma w_R(p)(1-r^*) - w_R(p)r^*(p+1), \\ \dot{V} = w_R(p)r^*V. \end{cases} \quad (\text{CRC-M}) $$

It can be seen that the study of the asymptotic behavior of system (S) applies to (CRC-M) as the latter is a particular case of the original one (S) with $r_0 = u^* = r^*$. We then maximize the final product $x^*$ in terms of the constant ribosomal concentration $r^* \in [r^-, r^+]$. The latter is given by the expression

$$ \max \lim_{r^* \to \infty} x(t). $$
We can see that the quantity

$$z = s + (p + 1)\mathcal{V} + \gamma \frac{1 - r^*}{r^*} \mathcal{V}$$

644 is constant. Thus,

$$\mathcal{V}^* + \gamma \frac{1 - r^*}{r^*} (\mathcal{V}^* - \mathcal{V}_0) = \Sigma$$

647 which, using (3.3), yields

$$x^* = \gamma \frac{1 - r^*}{r^*} (\mathcal{V}^* - \mathcal{V}_0).$$

650 Using the fact that $\mathcal{V}^*+x^* = \Sigma$ from (3.3), we see that $x^*$ is monotone decreasing w.r.t. $r^*$, and so the ribosomal concentration maximizing the infinite-time metabolite mass is $r^* = r^-$. This is what one would expect intuitively in an infinite-time horizon, as $r^* = r^-$ favors the production of metabolic proteins M, which catalyzes the synthesis of metabolites X without arresting the production of biomass (given by the case $r^* = 0$, which cannot be attained in trajectories starting in $\Gamma$). However, this kind of strategies might perform sub-optimally for the finite horizon case, as not having enough biomass can translate into a slow metabolite synthesis rate. Mathematically, this is represented through the presence of $\mathcal{V}$ in the dynamical equation of $x$. Similar to previous results [27], a first phase dedicated to bacterial growth can also foster the production of $X$, which depends directly on the concentration of bacteria in the bioreactor.

## 5.2. Finite-time problem.

### 5.2.1. Problem formulation.
In this section, we study the metabolite production objective in (S) for a time interval $[0, t_f]$, in which the final concentration of metabolite in the bioreactor $x(t_f)$ is maximized. While the biomass maximization objective $\mathcal{V}(t_f)$ was already studied in model (WTB-M) (representing a wild-type bacteria), it is likely that the presence of the heterologous pathway responsible for the production of $x$ might affect the results already obtained. Thus, the two objectives are compared in model (S) from a numerical perspective. Given a fixed final time $t_f > 0$, the OCP maximizing $c\mathcal{V}(t_f) + (1-c)x(t_f)$ (with $c = \{0, 1\}$ depending on the objective) with initial conditions (IC) writes

$$\begin{array}{ll}
\text{(MP-OCP)} & \left\{
\begin{array}{ll}
\text{maximize} & c\mathcal{V}(t_f) + (1-c)x(t_f), \\
\text{subject to} & \text{dynamics of (S),} \\
& \text{initial conditions (IC),} \\
& \text{and } u(\cdot) \in \mathcal{U},
\end{array}
\right. \\
672 & \\
673 &
\end{array}$$

It should be noted that the OCP is only valid for $c = 0$ (representing the metabolite production objective) and $c = 1$ (for the biomass maximization objective), and thus the intermediate values $c \in (0, 1)$ are not considered. One can easily see that, given the dynamics of the system, applying PMP would yield a Hamiltonian linear in the control for both values of $c$, which means that the solution of (MP-OCP) is similar to that of (BM-OCP), given by expression (4.1). However, (MP-OCP) has an additional
level of complexity in comparison with (BM-OCP), produced by the presence of $x$ in the model, as well as the term $-\gamma w_R(p)(1-r)$ in $\dot{p}$ responsible for the consumption of resources for metabolite production. Thus, it was not possible to perform a study of the OCP using PMP. A numerical analysis of these results is provided in the next section.

### 5.2.2. Numerical simulations.
The optimal trajectories were obtained following the same procedure as in the biomass maximization case. Figures 9, 10 and 11 are solutions of (MP-OCP) where the objective is the final-time product maximization $x(t_f)$ (obtained by fixing $c=0$) for different final times $t_f$ (40, 60 and 80, respectively) and with the same set of initial conditions. The metabolite synthesis rate is set to $\gamma = 0.5$. As expected, and similar to the results obtained for (BM-OCP), the optimal control takes the value $u=0$ for most of the interval, representing an allocation strategy that promotes the synthesis of proteins of the metabolic machinery $M$, consequently catalyzing the absorption of nutrients from the medium and the production of $x$. Solutions are characterized by a short $u=1$ bang arc at the beginning of the process, followed by a marginal singular arc before the final $u=0$ bang arc. The latter suggests that a valid sub-optimal approximation of the optimal control could be a simple bang-bang (1-0) control law. In that case, the only degree of freedom would be the switching time between bang arcs, that can be easily computed through numerical optimization methods. Additionally, and as it can be seen across Figures 9, 10 and 11, these results do not depend on the final time $t_f$: the final bang $u=0$ of the optimal control is always predominant in the control strategy, and becomes larger as $t_f$ increases. The finite-time numerical results are consistent with the results obtained for the infinite-time case in Section 5.1.1, in which the ribosomal sector of the cell $r$ should be minimized to maximize the production of $x$. Figure 12 shows an optimal trajectory solution of (MP-OCP) with cost function $V(t_f)$. In this case, the allocation strategy is described by an initial bang $u=1$ followed by a singular arc that takes up most of the optimal solution, with values near to an intermediate strategy $u=0.5$; and a short bang $u=0$ at the end. Such strategy leads to a bacterial composition much more balanced between ribosomal and enzymatic proteins, in opposition to the metabolite production case (with $c=0$ for (MP-OCP)), where most of the bacterial proteins were dedicated to the metabolic machinery. The latter behavior illustrates a natural trade-off between two opposed strategies: maximizing the number of ribosomes to prioritize the synthesis of macromolecules over the production of $x$ and, at the same time, maximizing the enzymatic activity in order to consume the substrate in the medium as fast as possible, towards (Full depletion).

Under the hypothesis that the mechanisms behind the allocation of cellular resources in bacteria have been optimized to outgrow competitors, it is coherent to think that a genetically modified bacteria (e.g. able to synthesize the metabolite X) would also maximize biomass. Thus, these internal mechanisms would produce a cellular composition profile similar to the one shown in Figure 12. Indeed, this is expected to happen even for artificially engineered specimens, as, in microorganisms, natural selection occurs very rapidly (even for time windows in the order of hours) due to the strong genetic variability of bacteria and their astonishingly high doubling rate. Then, interfering with this strategy so as to obtain allocations maximizing the production of metabolites in the bioreactor (as the ones shown in Figures 9, 10 and 11) can be accomplished by externally shutting down the production of ribosomes at a certain time instant, which can be triggered by well-known biotechnological control techniques such as growth arrest [11].
<figure><img src="image_10.png" /><figcaption>FIG. 9. Solution of (MP-OCP) for the metabolite maximization case $x(t_f)$, with $s_0 = 0.1$, $p_0 = 0.001$, $r_0 = 0.1$, $V_0 = 0.003$ and $t_f = 60$. The final product concentration $x(t_f)$ is at 35% of the total mass in the bioreactor $V_e\Sigma/\beta$, while the final volume $V(t_f)$ is only at 16%.</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>FIG. 10. Solution of (MP-OCP) for the metabolite maximization case $x(t_f)$, with $s_0 = 0.1$, $p_0 = 0.001$, $r_0 = 0.1$, $V_0 = 0.003$ and $t_f = 60$. The final product concentration $x(t_f)$ is at 62% of the total mass in the bioreactor $V_e\Sigma/\beta$, while the final volume $V(t_f)$ is only at 23%.</figcaption></figure>

**6. Discussion.** This paper presented a mathematical study of bacterial resource allocation in batch processing, and its applications to biomass and metabolite production. A dynamical model considering the production of a value-added chemical compound is proposed, and a study of the asymptotic behavior of the system based on mass conservation laws shows that, under all possible resource allocation strategies, all the substrate in the medium is consumed. Then, the particular case of a wild-type bacteria with no metabolite production is analyzed, showing that the optimal allocation propitious for biomass maximization—and thus, competitors outgrowing—is accomplished through a rather low value of the optimal control $u$, which yields a very high $m/r$ ratio (i.e. the ratio of enzymatic to ribosomal mass fractions) in the cell
<figure><img src="image_12.png" /><figcaption>FIG. 11. Solution of (MP-OCP) for the metabolite maximization case $x(t_f)$, with $s_0 = 0.1$, $p_0 = 0.001$, $r_0 = 0.1$, $V_0 = 0.003$ and $t_f = 60$. The final product concentration $x(t_f)$ is at 74% of the total mass in the bioreactor $V_e\Sigma/\beta$, while the final volume $V(t_f)$ is only at 19%.</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>FIG. 12. Solution of (MP-OCP) for the biomass maximization case $V(t_f)$, with $s_0 = 0.1$, $p_0 = 0.001$, $r_0 = 0.1$, $V_0 = 0.003$ and $t_f = 60$. The final product concentration $x(t_f)$ is at 27% of the total mass in the bioreactor $V_e\Sigma/\beta$, while the final volume $V(t_f)$ is at 56%.</figcaption></figure>

throughout the bioprocess. Paradoxically, for the metabolite production case, this kind of strategies would rather maximize the production of $x$, while for maximizing the biomass under the presence of the heterologous pathway, a more balanced cellular composition is required. Overall, results show that optimal allocation can be accomplished through a two-phase control: a first phase of pure bacterial growth dedicated to produce as many ribosomal proteins as possible; followed by a production phase where the remainder of the feedstock is used to synthesize the compound of interest. While the first phase matches the natural bacterial behavior, the second one requires human intervention to arrest the production of ribosomes (by externally
setting $u = 0$). The obtained two-phased external control profile is in close agreement with well-known metabolic engineering techniques used in microbial cell factory [27].

The analysis presented in this paper raises interesting questions both from math- ematical and biological points of view. For instance, it would be worthwhile to further study the potential presence (or absence) of the turnpike phenomenon in the optimal control solutions. Depending on the complexity of the OCP, it is often possible to obtain an analytical proof of the exponential convergence of the singular arc to the solution of the static OCP [25], and to find an explicit link between the length of the singular arc and the duration of the bioprocess. From a biological perspective, includ- ing additional substrates could broaden the approach to represent other well-studied phenomena observed in bacterial growth. For example, under the presence of differ- ent nutrients in the culture, bacteria tend to favor (i.e. consume first) those that are easier to metabolize, which is a phenomenon known as diauxic growth. This behavior has been extensively studied from an optimal control viewpoint, but without taking into consideration cellular composition. Formulating more comprehensive dynamical models and studying their associated OCPs can be instrumental in understanding natural allocation strategies in microbial growth, and in engineering synthetic control schemes targeting industrial objectives.

**Acknowledgments.** The authors thank Hidde de Jong (Microcosme team, Inria Grenoble – Rhône-Alpes) for many discussions on the biological context. We also appreciate the help of Sacha Psalmon and Baptiste Schall, from Polytech Nice Sophia, for the numerical simulations and the production of the online example in the control toolbox gallery.<fnref n="5" />

REFERENCES

[1] A. A. AGRACHEV AND Y. SACHKOV, *Control Theory from the Geometric Viewpoint*, vol. 87, Springer Science & Business Media, 2013.

[2] M. BASAN, M. ZHU, X. DAI, M. WARREN, D. SÉVIN, Y.-P. WANG, AND T. HWA, *Inflating bacterial cells by increased protein synthesis*, Molecular systems biology, 11 (2015), p. 836.

[3] V. BORISOV, *Fuller's phenomenon*, Journal of Mathematical Sciences, 100 (2000), pp. 2311–2354.

[4] J.-B. CAILLAU, W. DJEMA, J.-L. GOUZÉ, S. MASLOVSKAYA, AND J.-B. POMET, *Turnpike property in optimal microbial metabolite production*, J. Optim. Theory. Appl., 194 (2022), pp. 375–407.

[5] E. CINQUEMANI, F. MAIRET, I. YEGOROV, H. DE JONG, AND J.-L. GOUZÉ, *Optimal control of bacterial growth for metabolite production: The role of timing and costs of control*, in 2019 18th European Control Conference (ECC), IEEE, 2019, pp. 2657–2662.

[6] H. DE JONG, J. GEISELMANN, AND D. ROPERS, *Resource reallocation in bacteria by reengineering the gene expression machinery*, Trends in microbiology, 25 (2017), pp. 480–493.

[7] N. GIORDANO, F. MAIRET, J.-L. GOUZÉ, J. GEISELMANN, AND H. DE JONG, *Dynamical allocation of cellular resources as an optimal control problem: novel insights into microbial growth strategies*, PLoS computational biology, 12 (2016), p. e1004802.

[8] R. HEINRICH AND S. SCHUSTER, *The regulation of cellular systems*, Springer Science & Business Media, 2012.

[9] S. HUI, J. M. SILVERMAN, S. S. CHEN, D. W. ERICKSON, M. BASAN, J. WANG, T. HWA, AND J. R. WILLIAMSON, *Quantitative proteomic analysis reveals a simple strategy of global resource allocation in bacteria*, Molecular systems biology, 11 (2015), p. 784.

[10] L. HUO, J. J. HUG, C. FU, X. BIAN, Y. ZHANG, AND R. MÜLLER, *Heterologous expression of bacterial natural product biosynthetic pathways*, Natural Product Reports, 36 (2019), pp. 1412–1436.

<footnote n="5">⁵ct.gitlabpages.inria.fr/Gallery/substrate/depletion.html</footnote>
[11] J. IZARD, C. D. G. BALDERAS, D. ROPERS, S. LACOUR, X. SONG, Y. YANG, A. B. LINDNER, J. GEISELMANN, AND H. DE JONG, *A synthetic growth switch based on controlled expression of RNA polymerase*, Molecular systems biology, 11 (2015), p. 840.

[12] G. JEANNE, A. GOELZER, S. TEBBANI, D. DUMUR, AND V. FROMION, *Dynamical resource allocation models for bioreactor optimization*, IFAC-PapersOnLine, 51 (2018), pp. 20–23.

[13] F. MAIRET, J.-L. GOUZÉ, AND H. DE JONG, *Optimal proteome allocation and the temperature dependence of microbial growth laws*, NPJ systems biology and applications, 7 (2021), pp. 1–11.

[14] D. MOLENAAR, R. VAN BERLO, D. DE RIDDER, AND B. TEUSINK, *Shifts in growth strategies reflect tradeoffs in cellular economics*, Molecular systems biology, 5 (2009), p. 323.

[15] L. S. PONTRYAGIN, *Mathematical theory of optimal processes*, Routledge, 2018.

[16] M. SCOTT, C. W. GUNDERSON, E. M. MATEESCU, Z. ZHANG, AND T. HWA, *Interdependence of cell growth and gene expression: origins and consequences*, Science, 330 (2010), pp. 1099–1102.

[17] M. SCOTT, S. KLUMPP, E. MATEESCU, AND T. HWA, *Emergence of robust growth laws from optimal regulation of ribosome synthesis*, Molecular Systems Biology, 10 (2014), p. 747.

[18] I. S. TEAM COMMANDS, *Bocop: an open source toolbox for optimal control*. http://bocop.org, 2017.

[19] E. TRÉLAT AND E. ZUAZUA, *The turnpike property in finite-dimensional nonlinear optimal control*, Journal of Differential Equations, 258 (2015), pp. 81–114.

[20] A. Y. WEISSE, D. A. OYARZÚN, V. DANOS, AND P. S. SWAIN, *Mechanistic links between cellular trade-offs, gene expression, and growth*, Proceedings of the National Academy of Sciences, 112 (2015), pp. E1038–E1047.

[21] A. G. YABO, *Optimal resource allocation in bacterial growth: theoretical study and applications to metabolite production*, PhD thesis, Université Côte d’Azur, 2021, https://theses.hal.science/tel-03636842.

[22] A. G. YABO, J.-B. CAILLAU, AND J.-L. GOUZÉ, *Singular regimes for the maximization of metabolite production*, in 2019 IEEE 58th Conference on Decision and Control (CDC), IEEE, 2019, pp. 31–36.

[23] A. G. YABO, J.-B. CAILLAU, AND J.-L. GOUZÉ, *Optimal bacterial resource allocation: metabolite production in continuous bioreactors*, Mathematical Biosciences and Engineering, 17 (2020), pp. 7074–7100.

[24] A. G. YABO, J.-B. CAILLAU, AND J.-L. GOUZÉ, *Optimal allocation of bacterial resources in fed-batch reactors*, in 2022 European Control Conference (ECC), IEEE, 2022, pp. 1466–1471.

[25] A. G. YABO, J.-B. CAILLAU, J.-L. GOUZÉ, H. DE JONG, AND F. MAIRET, *Dynamical analysis and optimization of a generalized resource allocation model of microbial growth*, SIAM Journal on Applied Dynamical Systems, 21 (2022), pp. 137–165.

[26] A. G. YABO AND J.-L. GOUZÉ, *Optimizing bacterial resource allocation: metabolite production in continuous bioreactors*, IFAC-PapersOnLine, 53 (2020), pp. 16753–16758.

[27] I. YEGOROV, F. MAIRET, H. DE JONG, AND J.-L. GOUZÉ, *Optimal control of bacterial growth for the maximization of metabolite production*, Journal of mathematical biology, 78 (2019), pp. 985–1032.