Oscillatory behaviour control in a continuous culture under double-substrate limitation

ABSTRACT The paper discusses the possibility of using a mixture of two growth limiting substrates to induce or eliminate self-sustained oscillations in a continuous culture process. The proportion of both substrates in the mixture is treated as a new control variable. The presented approach is based on the assumption that the oscillatory behaviour occurs for selected substrates in some range of dilution rates. Because a double-substrate limitation may occur, the analysis is performed for two fundamental substrate utilization patterns: simultaneous consumption and diauxic growth. By using model simulations and bifurcation analysis, we show that an appropriate proportion of two substrates in the mixture allows for the control of the oscillatory behaviour.


Introduction
Recent years have witnessed an increasing interest in the theoretical and experimental analysis of continuous fermentation processes. Among many works on this subject, a large group of papers can be distinguished, where the authors try to explain the reasons for the occurrence of self-sustained oscillations (SSO) of biomass concentration in a batch [28] or continuous culture for constant dilution rate D and inlet substrate concentration S in [7,18,50,55]. However, in spite of many experimental studies, these reasons are still not known precisely [45].
In the literature, the processes of continuous ethanol fermentation with yeast Saccharomyces cerevisiae [5] or bacteria Zymomonas mobilis [6,40,57] limited for growth with glucose, are often discussed. The two most frequently mentioned causes of SSO are: a time delay in the response of cells to changes in the environment, which is caused by the adaptation of cells to new conditions [26,37], or the interaction between linked intracellular reactions, where the most commonly quoted example is the synchrony of cell division [13,19].
Apart from the descriptions of the mechanisms of SSO, the influence of oscillatory behaviour on the process performance is also studied. For example, it is possible to obtain greater values of biomass concentrations or biomass productivity in the bioreactor operating in the range of SSO in comparison to the results obtained at steady-states [8,42,51]. However, the operation of bioreactor in the range of SSO is usually avoided in practice. Especially, it concerns industrial-scale continuous bioreactors for which the oscillatory mode of operation can lead to high loss of residual substrates [16,57]. Hence, it is necessary to investigate the possibilities of induction or elimination of the SSO of biomass concentration in continuous fermentation processes. This can be achieved by changing the operating conditions (e.g. dilution rate, pH, temperature, dissolved oxygen level), as shown experimentally in [34,44]. However, in many cases, these parameters have to be kept constant, since they determine the optimal conditions for the growth of microorganisms. As a result, most bioreactors operate at steady-state conditions for small dilution rates [8], which may not be the optimal solution.
Another possibility is to use an additional substrate without changing the other process parameters. Its inlet flow rate or concentration can be treated as the additional control variable. In the literature, one can find examples of experiments that use additional substrates or secreted products to eliminate (or induce) the SSO in continuous fermentation processes [34,44]. For instance, in the experiments performed by Parulekar and co-workers [44], the presence of SSO was detected for a wide range of dilution rates D and for sufficiently high levels of dissolved oxygen. In these experiments, small amounts of ethanol were added to the oscillating culture causing an immediate arrest of oscillations that reappeared after consumption of the ethanol. In turn, an additional portion of glucose increased the amplitude of the oscillations which returned to the initial value after utilization or washout of the additional glucose. The effect of additional portions of the product (ethanol) for the existence of SSO was also investigated by Martegani and co-workers [34]. The conducted experiments showed that the injection of ethanol to the oscillating culture did not cause any changes when the cells were already producing ethanol, while a marked phase shift in the oscillations was obtained in the absence of ethanol in the culture. On the other hand, the addition of ethanol to the non-oscillating culture triggered new SSO. Furthermore, based on the numerical investigations, Sridhar [57] considered the case of adding small amounts of substrate, biomass, or product, which caused the disappearance of the oscillations in the system. The obtained results confirmed, to some extent, the possibility of the SSO control in the presence of two limiting substrates (e.g. ethanol and glucose). However, these experiments were carried out for finite portions of substrates or products and the case of continuous feeding of substrates was not taken into account. A similar numerical analysis was also performed in [40,41] by studying the influence of dilution rate and additional amounts of cells or ethanol to the feed stream on the elimination of SSO. In turn, Beuse and co-workers [13] experimentally investigated the influence of the type of carbon source on the mode of oscillations for a fixed dilution rate. Their results have shown that the change of carbon source (in their case galactose was replaced by glucose) prolonged the period of oscillations, but the case of using a mixture of these substrates was not investigated. The problem of using an additional substrate was also discussed in [52] for the classical chemostat model with variable yield coefficient and in [53] for the continuous fermentation process with a product inhibition model, where the mixture of two substrates was continuously fed into the reactor vessel. In both cases, the numerical results have shown that an appropriate proportion of both substrates in the mixture allowed for induction or elimination of the oscillatory behaviour. However, the results presented in [52] were obtained based on the assumption that both substrates were utilized simultaneously, while in fermentation processes, a diauxic growth pattern is frequently observed as well [20]. In turn, the results presented in [53] were obtained assuming that the growth of microorganisms is inhibited by high product (ethanol) concentrations and only the diauxic growth pattern was taken into account.
In comparison to the previous works, we investigate a more general case of the continuous culture system described by an unstructured model with constant and variable yield coefficients and including two fundamental substrate utilization patterns. The novelty of the proposed approach is the use of a mixture of two different substrates that can be utilized either simultaneously or sequentially according to the diauxic growth pattern. One of the main goals of the paper is to evaluate the usefulness of a new control variable (the proportion of both substrates in the mixture) in order to induce or eliminate the oscillatory behaviour. To reduce costs of the approach, the additional substrate will be fed into the bioreactor at the expense of the original substrate. The sum of both inlet flow rates is always kept constant and only their proportion is varied.
The paper is organized as follows: a detailed description of the mathematical model including single and double substrate limitation, and the idea of SSO control are presented in the next two sections. Section 4 contains the main results for two substrate utilization patterns and Section 5 concludes the paper.

An unstructured mathematical model for a single limiting substrate
The classical chemostat model with variable yield coefficient has the following form [22,54]: where: μ m -maximum specific growth rate [1/h]; K s -half saturation constant [g/L]. In order to simplify the discussion, an additional equation for the product concentration has not been taken into account. Moreover, the model (1)-(3) does not include the death of microorganisms; thus we omit the term -k d X (k d -death rate constant) in equation (2). In our analysis, the inclusion of this term does not change the basic behaviour of the system.
It is well known that the system (1)-(3) with the constant yield coefficient (Y = const) does not exhibit the oscillatory behaviour [2,3,58]. However, if the yield coefficient is dependent on the substrate concentration S: where α [g/g], β [L/g] are constant and positive coefficients, then there exists a stable limit cycle for a range of dilution rates D in the system (1)-(4) [1][2][3][4]42]. This is also consistent with the experimental results presented in [12,13]. Hence, two main cases are distinguished: Case a) the bioreactor is fed with a substrate for which there are no oscillations for any value of dilution rate D, and the yield coefficient is assumed to be constant. Accordingly, the substrate is referred to as the substrate not generating SSO (S ngo ). Case b) the bioreactor is fed with a substrate for which there may exist SSO for some range of dilution rates D, and the yield coefficient is given by (4). In this case, the substrate is referred to as the substrate generating SSO (S go ). Figure 1 presents the behaviour of the system (1)-(3) with a constant or variable yield coefficient for parameters values taken from [8] and shown in Table 1. The black dots represent the range of SSO and the grey shaded area represents maximum and minimum values of biomass concentration.
Of course, as it was mentioned in the previous section, it should be emphasized that the assumption that the yield coefficient is given by (4) does not explain the occurrence   of the SSO of biomass concentration, but provides the desired dynamical properties of the bioreactor system (1)-(3). Therefore, simulation results of the chemostat with variable yield coefficient (4) may not agree satisfactorily with experimental data in some cases [27]. On the other hand, a wide class of continuous culture processes can be represented by the model (1)-(4). The assumed dependency (4) of the yield coefficient Y on the substrate concentration S is commonly encountered in the literature [42,46,62] and results from the linearized form of the Pirt equation [24,47] that is used to estimate the maintenance energy requirements [38,47,60]. In other words, it is assumed that some fraction of the substrate is used to support basic life functions of microorganisms. This implies that the yield coefficient is dependent on the specific growth rate of microorganisms, which was also confirmed in experimental studies [15].

The idea of using a mixture of two growth limiting substrates
Assuming that the change of substrate (carbon source) may significantly affect the dynamical properties of the process, it will be shown that the induction or elimination of SSO of the biomass concentration is possible by setting an appropriate proportion of the individual substrates (S ngo and S go ) in the mixture, fed at a constant flow rate Q into the bioreactor. The presentation and discussion of this possibility is the main thesis of this paper, since the proportion of substrates can be treated as the additional control variable.
A general scheme of the bioreactor fed with the mixture of two different substrates of inlet concentrations S 1in and S 2in is presented in Figure 2. The idea of modelling the use of two different feed vessels, each with a different growth-limiting substrate and input rate, was first introduced in [17].
The proportion r∈[0,1] of both substrates in the inlet stream is determined by means of the three-way valve and the individual flow rates for the first substrate (S 1in ) and for the second substrate (S 2in ) are equal to (1-r)Q and rQ, respectively. Thus, the flow rate for the mixture of both substrates is always equal to (1-r)Q + rQ = Q.

Modelling of the continuous culture with two growth limiting substrates
The bioreactor is fed with the mixture of two different substrates ( Figure 2) and if both substrates concentrations are below their saturation levels, then a double-substrate limitation may occur. Hence, the general form of the mathematical model is as follows:  (7), since it does not change the basic behaviour of the system. According to the remarks presented in the preceding section, it is assumed that the yield coefficient Y 1 is linearly dependent on the substrate concentration S 1 and the yield coefficient Y 2 for the additional substrate, is constant. As far as the form of the specific growth rate μ (8) is concerned, there are many different models for μ given in the literature. In general, they can be classified into two groups: non-interactive and interactive doublesubstrate limited growth models [25,32,49,59,63]. Thus, finding a model is simply reduced to finding a function for the specific growth rate μ (8). In the case of non-interactive models, only one of the two substrates is limiting at a time [25,63]. In this case, some authors assumed that the specific growth rate is equal to the lowest growth rate predicted from the separate single-substrate models [17]. In turn, for the interactive models (where the growth rate of microorganisms is limited by two substrates at a time), the function for μ is more complex and it may have one of the following forms: product [56], non-weighted sum [11,32,61] or weighted sum of individual specific growth rates for each of the substrates [14,30,33]. However, it should be emphasized that most of these models were tested mainly for batch cultures, and only few of them have been verified in laboratory experiments for the continuous mode of operation [31,56,61]. This is the consequence of the lack of experimental data on the mixture of substrates [29], especially for small dilution rates D [31]. A reader interested in a more detailed review of models for the mixture of two or more substrates should refer to [29,31,43].
On the other hand, it should be mentioned that the growth of microorganisms in the presence of two sources of substrates is a quite complex process, as the available substrates can be utilized according to various consumption patterns. The presence of several limiting substrates can also affect the dynamical properties of continuous culture processes [9,10,17,25,36,52,53]. Egli [20] presented a thorough study of the possible consumption patterns in the presence of several available substrates. Among these patterns, the most frequently observed are: simultaneous utilization of both substrates in the whole range of dilution rates (until washout occurs) and the diauxic growth [20,27,48]. The diauxic phenomenon in the continuous culture with two growth limiting substrates: the more-preferred ( ) and the less-preferred substrate ( ) [21]; S 1in and S 2in are the inlet substrate concentrations for the more-preferred and the less-preferred substrate, respectively. Figure 3 presents the theoretical steady-state diagram for biomass and residual substrate concentrations, in the case of diauxic growth pattern, which is observed when both substrates are utilized for relatively small dilution rates D, i.e. for D ≤ D trans (where D trans is a transition dilution rate). In turn, for D > D trans the only utilized substrate is a substrate on which microorganisms can grow faster [20,21]. The observed behaviour is a consequence of the fact that microorganisms maximize their growth rate. In the continuous culture at steady-state conditions the specific growth rate of microorganisms is equal to the dilution rate D. Hence, for a sufficiently large dilution rate (D > D trans ), the microorganisms utilize the substrate that is able to support the enforced growth rate. The other substrate is washed out of the reactor vessel or its concentration is close to its inlet value. In the remainder of the paper, the substrate supporting the faster growth of microorganisms will be referred to as the more-preferred substrate, and the other one -the less-preferred substrate (Figure 3).
To study the effectiveness of the proposed idea of SSO control under various substrate consumption patterns, we assume that the specific growth rate μ is modelled by the generalized Monod equation proposed by Yoon and co-workers [61], which means that both substrates are perfectly substitutable [25,63]: where: a 1 = K s2 /K s1 , a 2 = K s1 /K s2 are dimensionless coefficients determining the inhibition effect of the individual substrates on the growth of microorganisms. Thus, assuming that S 1in , S 1 are the concentrations of the substrate generating SSO (S go ) and S 2in , S 2 are the concentrations of the substrate not generating SSO (S ngo ), the equations of the unstructured model are as follows: where: X 0 , S 10 , S 20 are initial conditions for biomass and substrates concentrations [g/L], respectively. Unlike other models discussed in the literature, the model (10)-(12) has the following important property. If the bioreactor is fed with the mixture of the same two substrates, i.e. K s1 = K s2 , μ m1 = μ m2 , a 1 = a 2 = 1, and the yield coefficients are constant Y 1 = Y 2 = const, the model (10)- (12) can be reduced to the classical chemostat system (1)-(2) with the specific growth rate given by the Monod relation (3), where the substrate concentration is equal to S: = S 1 + S 2 and the inlet substrate concentration is equal to S in : = (1−r)S 1in + rS 2in for any r∈[0,1].
Depending on the degree of opening of a three-way valve (Figure 2), two extreme cases can be distinguished. The former case is for r = 0, when the bioreactor is fed with the first substrate (S go ) of inlet concentration S 1in , which means that there are SSO for some range of dilution rates D. In the latter case, i.e. for r = 1, only the second substrate (S ngo ) of inlet concentration S 2in is fed to the bioreactor and there are no SSO for any value of D. However, depending on the proportion r∈[0,1], both steady-states and SSO (stable limit cycle) may occur. Although the mathematical model contains only three nonlinear equations, the analytical determination of all the equilibrium points is a complex task. It is possible to show that the system (10)- (12) can have multiple steady-states with hysteresis in a very narrow range of dilution rates. Therefore a major part of the analysis (including simulations and determination of bifurcation diagrams) will be performed in Matlab and XPPAut software [23] for parameter values shown in Tables 1 and 2. The parameter values in Table 1 have been chosen according to [8] and are referred to as the base parameters. In turn, Table 2 contains only maximum specific growth rates and half-saturation constants for the more-preferred and less-preferred substrates, and are referred to as the modified parameters. The parameters values have been chosen so as to reflect the diauxic growth pattern.
It is easy to notice that one of the equilibrium points is E = (S * 1 , S * 2 , X * ) = ((1 − r)S 1in , rS 2in , 0), which corresponds to the washout state. From the local stability analysis of this point, we determine a critical value of dilution rate Dc. The Jacobian matrix at the washout equilibrium point is: and eigenvalues of (13) are: Then E = ((1−r)S 1in , rS 2in , 0) is locally asymptotically stable if D > D c , where D c is the critical value of dilution rate: It is possible to show that the washout equilibrium point E is globally asymptotically stable in Ω, if the dilution rate D is sufficiently large. By using (9), it can be shown that μ 1 (S 1 ,S 2 ) + μ 2 (S 1 ,S 2 ) < max{μ m1 , μ m2 } for any S 1 ≥ 0 and S 2 ≥ 0. Let μ m1 ≥ μ m2 , then from (9) we have: Now, multiplying inequality (16) by K s1 K s2 + S 1 K s2 + S 2 K s1 > 0, we get: Hence, if D ≥ max{μ m1 , μ m2 }, then X(t) tends to 0 as t tends to infinity for any initial conditions (S 10 ,S 20 ,X 0 )∈ Ω. In fact, the condition D ≥ D c does not always guarantee that the inequality μ 1 (S 1 ,S 2 ) + μ 2 (S 1 ,S 2 ) < D c is true for any S 1 ≥ 0 and S 2 ≥ 0. As has been mentioned earlier, the system (10) ) and the steady-state biomass concentration is dependent on the initial conditions (S 10 ,S 20 ,X 0 )∈ Ω. The equilibrium points for no washout conditions can be found from the set of following equations: The above equations are obtained by setting all the derivatives in (10)-(12) equal to zero. Then, (19) results from (12) for μ(S 1 ,S 2 ) = D and (17) results from (10) and (11) knowing that μ 1 (S 1 ,S 2 ) + μ 2 (S 1 ,S 2 ) = D. Finally, (18) can be easily obtained from (10) and (11) noticing that μ 1 (S 1 ,S 2 )·S 2 K s1 μ m2 = μ 2 (S 1 ,S 2 )·S 1 K s2 μ m1 .
Further investigation of the system (10)-(12) is focused on studying the influence of the dilution rate D and the proportion r on the existence of SSO. In the case of the occurrence of SSO, the average biomass concentration is plotted on the graphs. First, we investigate the case of simultaneous utilization of both substrates and, in the next subsection, the diauxic growth pattern.

Simultaneous utilization of substrates
In this case, the model parameter values were taken from Table 1. For simplicity, but without loss of generality, we assume the same parameter values for the maximum specific growth rates μm 1 = μm 2 and for the half-saturation constants Ks 1 = Ks 2 . The simulation results in Figure 4 show that both substrates are consumed in the whole range of dilution rates until washout occurs.
An increase of the residual substrate concentration S go for D > 0.15 [1/h] and r = 0.25 results from changes in the yield coefficient Y 1 (S 1 ) which in turn changes the uptake rate for the substrate S go .
The main results are shown in Figure 5. The black thick lines represent the extreme cases for r = 0 and r = 1, and the thin black lines represent the maximum and minimum values of the biomass concentrations in the presence of SSO for a chosen value of r. By feeding   the bioreactor with additional substrate S ngo , the range of dilution rates D in the presence of SSO is becoming narrower and, eventually, disappears for sufficiently high values of r ( Figure 6).
As shown in experimental studies by Bruce and co-workers [16], the oscillatory behaviour was observed only for sufficiently high inlet substrate concentrations. This may explain the results shown in Figure 6, since for sufficiently low inlet concentrations of the S go substrate (higher values of r) the oscillations of the biomass concentration disappeared completely. As a result, three different regions can be distinguished: the region of SSO, the region of steady-states and the region of washout, which, irrespective of the proportion of both substrates, is always present for higher dilution rates ( Figure 6). Because in our case μ m1 = μ m2 = μ m , K s1 = K s2 = K s and S 1in = S 2in = S in , the critical dilution rate D c given by (15) is independent on the proportion r ( Figures 5  and 6): The simulation results have been complemented by the bifurcation diagrams obtained in XPPAut. Assuming that the dilution rate D is a bifurcation parameter, Figure 7 presents the steady-state biomass concentration for a few proportions r of substrates in the mixture. The thin continuous lines represent stable branches and the broken lines are unstable ones. The upper stable branch corresponds to an equilibrium point of nonzero biomass concentration. In some range of dilution rates, the system (10)-(12) becomes unstable and undergoes a Hopf bifurcation (Figure 7(a-c)).
It was also possible to find the regions of SSO (stable limit cycle), steady-states and washout on the parameter plane (D,r). The regions obtained in XPPAut agree with those presented in Figure 6.

Diauxic growth pattern
First, we assume that the substrate generating SSO (S go ) is the more-preferred one. Then, conversely, we assume that the substrate not generating SSO (S ngo ) is more-preferred by the microorganisms.
For these two cases, the maximum specific growth rates and the half-saturation constants have been chosen according to Table 2 and the other parameter values according to  Table 1. Depending on which substrate is the more-preferred one, Figures 8 and 9 show that for small dilution rates both substrate are consumed simultaneously and for higher dilution rates, the more-preferred substrate is mostly utilized by microorganisms.
The exceptions are the results obtained for r = 0.25 (Figure 8), which is a consequence of strong variations in the yield coefficient Y 1 (S 1 ). In this case, one can observe a decrease in the uptake of substrate generating oscillations (S go ) for D > 0.18 [1/h].
The influence of D and r on the existence of SSO has been presented in Figure 10 and the characteristic regions of the SSO, steady-states and washout -in Figure 11. As can be   clearly seen, irrespective of the substrate utilization pattern, an appropriate proportion of both substrates in the mixture allows for control of the oscillatory behaviour.
However, if the S go substrate is the more-preferred one, the region of SSO is larger and shifted toward higher dilution rates ( Figure 11). In the opposite case, the region of SSO is smaller and shifted toward lower dilution rates ( Figure 11). This can be explained by the fact that for a transition dilution rate the microorganisms switch their substrate utilization pattern (Figures 8 and 9) and the only utilized substrate is a substrate on which microorganisms can grow faster [20].
All the obtained results have been checked and complemented by the bifurcation diagrams obtained in XPPAut for two chosen values of r (Figures 12 and 13).   The behaviour of system for a single limiting substrate (r = 0 or r = 1) has already been shown in Figures 10 and 11. In these two cases, due to the differences in the maximum specific growth rates and the half-saturation constants, the washout phenomenon may occur for smaller dilution rates ( Figure 11) and the critical dilution rate D c can be calculated by using (15).
Similar results were obtained in [53] for the mixture of two substrates that were utilized according to the diauxic growth pattern with the specific growth rate μ given by (9). The main assumption was that the oscillatory behaviour resulted from the product inhibition and delayed response of microorganisms to changes in the product (ethanol) concentration. Despite this fact, it was possible to induce or eliminate the oscillatory behaviour by changing the proportion of both substrates in the mixture. At the same time, similar characteristic regions, as those shown in Figure 11, were also obtained [53].

Concluding remarks
The presented results were obtained for the continuous culture of microorganisms with constant and variable yield coefficients under double substrate limitation. Irrespective of the substrate utilization pattern and depending on the proportion r of two substrates in the mixture, the continuous culture system exhibited steady-state or oscillatory behaviour ( Figure 14). Qualitatively the same results were obtained in [53] for the continuous fermentation process described by a different mathematical model including the product inhibition and delayed response of microorganisms to changes in the product concentration.
The system (10)-(12) with constant yield coefficients was studied in [9,10] and the presented analysis did not reveal the oscillatory behaviour in the system. The assumption that one of the yield coefficients is variable can greatly change the dynamical behaviour of the chemostat system. As shown in [54], the chemostat with two species competing for a single substrate and with constant yield coefficients does not exhibit the oscillatory behaviour. Moreover, the competitive exclusion principle holds, which means that at most one population of microorganisms can survive and the losing population is washed out of the reactor vessel. Arino et al. [4] have studied the chemostat system with two species competing for a single substrate with specific growth rates given by Monod equation (3), and with a variable and constant yield coefficients. This means that one of the uptake functions was decreasing at high substrate concentrations and the other one was increasing. The numerical bifurcation analysis has shown that by introducing microbial population with a variable yield coefficient resulted in the occurrence of stable limit cycle for some range of parameter values. In effect, two competing populations could coexist.
In our case, the introduction of substrate for which the yield coefficient is variable and the specific growth rate function is given by the generalized Monod equation (8), also led to the occurrence of stable limit cycles. In both cases, i.e. for chemostat systems with a mixture of two substrates or with two competing species, a common assumption was that the system with a single microbial species and single growth limiting substrate exhibited oscillatory behaviour. Arino et al. [4] have shown that in the classical chemostat system (1)-(2) with a variable yield coefficient, a Hopf bifurcation can occur if the specific growth rate function μ(S) is increasing and the uptake function u(S): = μ(S)/Y(S) is decreasing at high substrate concentrations. These hold for the specific growth rate given by Monod equation (3) and the yield coefficient given by (4).
Hence, the additional substrate fed at a constant flow rate into the bioreactor, introduces an additional degree of freedom to control SSO, which means that the key process parameters (e.g. temperature, dilution rate D or inlet substrate concentrations) can be kept constant. As a result, the proportion r can be treated as an alternative or additional control variable in designing of a feedback controller that is robust to disturbances or uncertainties in model parameters, as shown for example in [35].