Computer Simulation of Bioprocess

Bioprocess optimization is important in order to make the bioproduction process more efficient and economic. The conventional optimization methods are costly and less efficient. On the other hand, modeling and computer simulation can reveal the mecha- nisms behind the phenomenon to some extent, to assist the deep analysis and efficient optimization of bioprocesses. In this chapter, modeling and computer simulation of microbial growth and metabolism kinetics, bioreactor dynamics, bioreactor feedback control will be made to show the application methods and the usefulness of modeling and computer simulation methods in optimization of the bioprocess technology.


Introduction
Bioindustry is important in utilization of reproducible resources, developments of environmental friendly production processes, and sustainable economy. In order to make the bioprocesses more efficient and economic, bioprocess optimization and automatic control are needed. The conventional optimization methods cost much labor, time, and money; on the other hand, modeling and computer simulation method can reveal the mechanisms behind the phenomenon to some extent, to assist the deep analysis and optimization of bioprocesses. The modeling and computer simulation method is much efficient and economic, and widely used in research and modern bioindustries.
Bioprocess efficiency depends on the cell capability, bioreactor performances, and the optimal control of the cultivation conditions. The metabolic network inside the cells involves thousands of enzymes, and the enzyme expression and activities are dynamically affected by the cultivation conditions. As a result, the cultivation condition affects the cell growth, metabolism, and product production in a sophisticated and nonlinear way. Control and maintain relatively optimal cultivation conditions through proper operation and control of the bioreactor are needed to improve the production efficiency of the bioprocess.
Bioprocess mathematical modeling involves the modeling of the dynamic changes of the metabolic rates and their distribution inside the cells with the changes of time and cultivation conditions, the modeling of the dynamic changes of the reaction rates and mass transfer rates as well as the cultivation conditions inside the bioreactor, and the modeling of the dynamics of the bioreactor control system etc., based on which optimizations of the bioreactor operation and control strategies can be made and the results can be predicted and evaluated by computer simulation. In this chapter, examples of modeling and computer simulation of microbial growth and metabolism kinetics, bioreactor dynamics, and the feedback control of the bioreactor are given to show the application methods and the usefulness of modeling and computer simulation methods in bioprocess technology.

Modeling of microbial cell growth and metabolism 2.1. Modeling of microbial cell growth
Cell growth is one of the most important variables to be investigated in bioprocess. The cell growth is usually described by the specific growth rate, μ, and the time course of cell concentration, X. The specific growth rate is defined by the increase in grams of cells (g) per gram dry cells (g) per hour (h), and can be modeled by Eq. (1) The specific growth rate is related with many process variables, like temperature (T), pH, dissolved oxygen (DO) concentration (C L ), substrate concentration (S), product concentration (P), X, and time (t). expressed by In real applications, only the key process variable(s) are included in Eq. (2) for simplification. Monod equation [1] using the substrate concentration as the single independent variable is shown by Eq. (3), as T, pH, and in many cases C L are controlled constant and can be neglected from the equation.
where μ m , is the maximum specific growth rate and k m , is the substrate affinity coefficient.
The typical cell growth curve is of "S" type, which has a lag growth phase and cannot be properly modeled by Monod equation as discussed later. At the initial cultivation stage, the cells need some time to adapt to the new environmental conditions for induction of some new enzymes needed for cell metabolism, etc., and the specific growth rate is zero or at a low value resulting in the lag growth phase. One way to model the lag growth phase is to separate the newly inoculated cells as active cells, X, and inactive cells, Y, and the time for Y to turn into X conforms to Pearson distribution expressed by Eqs. (4)-(7) [2].
where a 1 , a 2 , m 1 , and m 2 are the constants of the Pearson distribution. The other methods used to predict the lag growth phase simply defined the lag time in terms of cell growth, t L [2,3]. One way to deal with the lag growth phase is to relate it with the changes of μ defined by Eq. (8) [3].
After the lag growth phase, the specific growth rate increases gradually and the cells go into the exponential growth phase, which is expressed by Eqs. (9) and (10) where X 0 is the initial cell concentration. After the exponential growth phase, the specific growth rate decreases gradually to zero, because of nutrients limitation, accumulation of intracellular toxic intermediates, accumulation of inhibitors in the culture broth, etc., and the net cell growth tends to cease to enter into the stationary growth phase. In order to model the decreased cell growth rate and the stationary growth phase, the Logistic growth model was developed [4], in which μ decreases with the increase of cell concentration, X, and μ reaches zero when x reaches its maximum value, X m , shown by Eq. (11) From above analysis, it can be seen that the specific growth rate will start from zero or a low value in the lag growth phase, increases gradually and reaches the maximum value in the exponential growth phase, and then decreases gradually in the declined growth phase, which makes the time course of the specific growth rate the "bell" type curve and the time course of cell concentration the typical "S" type curve (Figure 1), which cannot be well fitted by the models discussed above. In order to simulate the "bell" type specific growth rate curve and the "S" type cell growth curve more accurately, the following model is developed shown by Eqs. (12) and (13) [5] μðtÞ ¼ μ m Á 1 1 þ e ÀkinðtÀtinÞ Á 1 1 þ e kdeðtÀtdeÞ ð12Þ dX dt ¼ μðtÞ Á X ð13Þ Figure 1. Graphically illustration of model parameters [5]. I is the lag growth phase; II is the increased growth phase; III is the exponential growth phase; IV is the decreased growth phase; V is the stationary growth phase. k in is the maximum increasing rate of μ; k de is the maximum decreasing rate of μ; t in is the time point when dμ/dt equals k in ; t in is the time point when dμ/dt equals k de ; and t L is the lag time.
where k in is the maximum increasing rate of μ; k de is the maximum decreasing rate of μ; t in is the time point when the increasing rate of μ equals k in ; t de is the time point when the decreasing rate of μ equals k de . All the parameters used in the model can be obtained graphically ( Figure 1).
One example of this model application is shown in Figure 2. In order to make wider application of above model, Eq. (12) can be combined with Monod model to develop Eqs. (12) and (13) into Eqs. (14) and (15) Figure 2. Simulation of cell growth of Trichoderma reesei [5].
Even if Monod model has some limitations, it is still the widely used growth model in real applications for the major reasons of simplification and the single independent variable of substrate concentration, which is the key process variable to be investigated in many fermentation processes. In cases of high density fermentation, or substrate or product inhibition, modifications of Monod model are needed. Contois model shown by Eq. (16) is an example for high density fermentation, in which modeling the cell concentration is included in the denominator of the specific growth rate equation to show the limitation effect of high cell concentration on the growth, to make the specific growth rate to be reciprocal to the cell concentration (μ ∝ X À1 ) at very low substrate concentration.
In some cases, the substrates which have inhibitory effect on cell growth, like ethanol or acetate, etc., are used. One example of the growth model under noncompetitive substrate inhibition with K I >> K m is shown by Eq. (17) One example for modeling product inhibition, like ethanol or lactic acid fermentation, is shown by Eq. (18) In case of dual substrates, the growth model in form of the sum or product of two Monod type terms is often used for the substitutable and nonsubstitutable substrates, respectively. For example, glucose and glycerol are substitutable substrates which can be modeled by Eq. (19), while glucose and oxygen are nonsubstitutable substrates which can be modeled by Eq. (20).
Above growth models are relatively simple, which are unstructured and unsegregated models, and are useful for practical applications. Structured and segregated growth models, which involve the intracellular structure or the nonhomogeneity of the cells, respectively, are generally sophisticated and contain uneasily measurable model parameters, and are usually used for theoretical purposes.

Modeling of microbial substrate uptake and product production
In many cases, the ratio of cell mass produced per substrate utilized is a constant and defined as the cell yield from the substrate, Y X/S , shown by Eq. (21) The minus sign in Eq. (21) is to ensure Y X/S to be positive as ΔS is negative. From Eq. (21), it can be seen that substrate consumption is proportional to the cell growth, so that substrate consumption can be simply modeled by Eq. (22) Further, the total substrate consumed can be considered of two parts, with one part for real cell growth and the other part for life maintenance to develop Eq. (22) into Eq. (23) where Y G is the maximum cell yield when μ tends to μ m ; m s the maintenance coefficient. From Eq. (23), Eq. (24) can be obtained showing the positive relationships between the specific growth rate, μ, and the cell yield, Y X/S .
From Eq. (24), it can be seen that in order to increase the cell yield, Y X/S , a high value of μ should be maintained.
The specific substrate consumption rate is defined by the consumption in grams of substrate (g) per gram dry cells (g) per hour (h), and can be modeled by Eq. (25) From Eqs. (23) and (25), Eq. (26) can be obtained Then, Eq. (22) or (23) can be expressed in a simple way by Eq. (27) For the metabolism of facultative anaerobes grown in oxygen limited condition, the cell yield varies greatly depending on the degree of oxygen limitation. Catabolism of 1 mole of glucose can produce 36 (or 38) mole ATP under aerobic condition or produce 2 mole ATP under anaerobic condition. The ATP-based cell yield, Y ATP , can be regarded a constant of 10 g drycell/mol ATP. So, the cell yield of Y X/S under anaerobic condition will only be 1/18 (or 1/19) of that under aerobic condition. The partition of the carbon source between aerobic pathway and anaerobic pathway in the catabolism is determined by oxygen supply or degree of oxygen limitation. Examples of cell growth under oxygen limited condition will be given in Sections 4.2 and 4.3.
In modeling of specific product production rate, Luedeking-Piret equation is most often used for its simplification and usefulness, which relates the specific product production rate to the growth related and nongrowth related parts by using α and β terms, respectively, as described by Eq. (28). The total product production is described by Eq. (29).

Modeling of bioreactor with different operation methods
Continuous stirred tank reactor (CSTR) is the most popular type of bioreactor, which can be operated in batch, fed-batch, and continuous modes. For batch culture, no substrate is fed into the bioreactor except air for aeration or acid or base for pH control, and no culture broth is taken out of the bioreactor during the fermentation process. For modeling of a typical batch culture, the specific rates of cell growth (μ), substrate uptake (q S ), and product production (q P ) introduced in Section 2, and the mass balance equations of Eqs. (30)-(32) can be used.
For fed-batch culture, substrate is fed into the bioreactor but no culture broth is taken out during the fermentation process, so that the liquid volume is increasing. For modeling fedbatch culture, V is variant and the mass balance equations of Eqs. (33)-(36) can be used. The specific rates of μ, q S , and q P introduced in Section 2 can be used in the modeling.
where F, the substrate feeding rate. Eqs. (33)-(35) can be transformed into Eqs. (37)-(39) In fed-batch culture, F can be continuous, for example, to be constant, linear increase, exponential increase with time, or uncontinuous, for example, operated in a repeated pulse-fed mode. Fed-batch culture has many advantages over batch culture. It has higher substrate conversion yield, extends production phase and can eliminate substrate inhibition or Crabtree effects, etc., and is widely used in industry.
Continuous culture is another kind of bioreactor operation method, with which method substrate is continuously fed into the bioreactor meanwhile the culture broth is continuously taken out of the bioreactor at the same rate so that the liquid volume remains unchanged. Continuous culture has the advantage of high production efficiency but the disadvantages of low substrate conversion yield, strain deterioration, and easy contamination, and is not often used in industry. As a result, examples of only batch and fed-batch cultures are investigated in next section.

Mathematical modeling
Bioethanol is produced by anaerobic fermentation using Saccharomyces cerevisiae, which can grow anaerobically through fermentative pathway (glycolysis) catabolizing 1 mole of glucose and producing 2 moles of ethanol and 2 moles of ATP. S. cerevisiae can also grow aerobically through tricarboxylic acid (TCA) cycle catabolizing 1 mole of glucose producing 6 moles of CO 2 and 38 moles of ATP. The cell yield from ATP, Y ATP , is relatively constant, which is about 10 g dry-cell mass/mole ATP. S. cerevisiae will grow much faster aerobically than anaerobically for the reason to have more ATP used for cell growth.
Fermentation period, which can be roughly divided into growth phase and production phase, is one major factor affecting the production cost. Fermentation period will be shortened if the cell growth phase is shortened. By employing an aerobic condition during the cell growth phase to fasten the cell growth and an anaerobic condition during the ethanol production phase, the fermentation period should be shortened while the ethanol production remained. The growth phase aerobic pulse stimulated ethanol fermentation and the normal anaerobic ethanol fermentation operated in batch mode are investigated and compared by modeling and simulation [6].
The specific glucose consumption rate (q S ) subject to substrate and product inhibition effects is modeled by Eq. (40). In Eq. (41), Q is the on-off switch between anaerobic (Q ¼ 1) and aerobic (Q ¼ 1) conditions (Q is not Q). Eq. (42) describes the ATP production from glucose under anaerobic or aerobic condition. The cell growth is based on the net ATP for cell synthesis shown by Eq. (43). Under aerobic condition, 6 moles of O 2 are needed for oxidizing 1 mole of glucose shown by Eq. (44). Ethanol is produced during the anaerobic production phase shown by Eq. (45).
where S and P are the glucose and product (ethanol) concentrations, respectively; P cri is the critical value of ethanol concentration for inhibition of glucose consumption; q S and q S.max are the specific glucose consumption rate and its maximum value, respectively; k s , ki s , and α are the constants; q ATP , q P , q O2 , and μ are the specific rates of ATP and ethanol productions, oxygen consumption, and the specific growth rate; m S.ATP is the ATP consumption constant for cell maintenance; Y ATP is the cell yield from ATP. The mass balance equations are shown by Eqs. (46)-(49) where OUR is oxygen uptake rate.

Simulation results
Simulation was made using above mathematical model (Figure 3). The references for the parameter values used in the model or in calculation of the parameter values used in the model are shown in Table 1. In aerobic condition, q S.mas of 1 g/g/h and m S of 0.1 mole/g/h, respectively, are used, which will be decreased and increased, respectively, compared with anaerobic condition.
The simulation results of the conventional anaerobic ethanol fermentation and the early growth phase aerobic pulse stimulated ethanol fermentation processes are shown in Figure 3.
In simulation of the early growth phase aerobic pulse stimulated ethanol fermentation process, Q was set to zero (Q ¼ 1) for the first 3 h. The q ATP and μ were abruptly decreased and the q S and q EtOH were abruptly increased with the shift of Q from 0 to 1 ( Figure 3B). The results showed that early stage aerobic pulse stimulated ethanol fermentation and had the advantage in shortening the fermentation period for more than 10 h compared with the conventional anaerobic ethanol fermentation.

Mathematical modeling
In this section, glucose feeding control based on dissolved oxygen (DO) will be investigated by using a fed-batch fermentation process using Escherichia coli, which is often used as the host for recombinant protein production. This control strategy, which relates glucose concentration with DO changes [10,11], is practical as DO sensor is widely used in fermentation technology. In the fermentation process, oxygen is continuously transferred into the liquid phase from the gas phase at a certain oxygen transfer rate (OTR) under aeration and agitation conditions; meanwhile, oxygen is continuously consumed by microbes at a certain oxygen uptake rate (OUR). After the cell reach high concentration, oxygen limitation occurs when OUR becomes larger than OTR so that DO decreases to nearly zero. On the other hand, E. coli catabolizes glucose aerobically through tricarboxylic acid (TCA) cycle, catabolizing one molecule of glucose into six molecules of CO 2 consuming six molecules of O 2 . When glucose is depleted, O 2 consumption stops (OUR ¼ 0) while OTR is positive so that DO rises suddenly. So, the sudden DO rise can be the indicator for glucose depletion and used as the signal for glucose feeding. After glucose feeding, glucose consumption and oxygen uptake resume and DO drops again. Then, the control system will monitor the next sudden DO rise for glucose feeding control, which strategy can maintain glucose concentration in low level and is called DO stat control strategy. This control strategy has been analyzed by modeling method [12]. The specific growth rate is modeled by Logistic equation shown by Eq. (50). The specific glucose consumption rate is modeled to include two parts, one for the net growth and the other for the maintenance shown by Eq. (51). The mole specific oxygen consumption rate is six times of the specific glucose consumption rate as shown by Eq. (52). The specific product production rate is modeled by Luedeking-Piret [Eq. (53)]. OUR and OTR are shown by Eqs. (54) and (55), respectively.
where q O2 and q P are the specific rates of O 2 consumption and product production, respectively; α, β, are the constants for Luedeking-Piret equation; k La , is the volume oxygen transfer rate; C L and C* are the dissolved oxygen concentration and its saturated value, respectively; M O2 and M Gluc are the molecular weights of O 2 and of glucose, respectively.
The mass balance equations for fed-batch culture can be made and transformed into Eqs. (56)-(60).

Simulation results
In the simulation, glucose pulse feeding was made when DO increased over 10% in order to avoid noise interruptions. In each pulse feeding, a dosage equivalent to 20 g/L increase in glucose concentration was fed. The initial glucose was depleted at about 75 h and the product concentration was a little over 6 g/L at that time. By glucose pulse feeding, the product concentration was more than doubled ( Figure 4). Glucose and DO concentrations go up and down in turn and fluctuate during the control period. By using this control strategy, glucose concentration can be maintained in an averaged low concentration, which is desired and helps to overcome the glucose effects and increase the product yields. In addition, the DO stat control strategy does not need the extra sensor and is easily applied.

Mathematical modeling
DO control is important in fermentation process. The level of DO can affect the metabolic flux distribution and the product yield and production efficiency. As oxygen has low solubility in water, DO control is a hard task for fermentation process. Compared with feedback control, DO feedforward-feedback (FF-FB) control has the advantage in dealing with the time-varying characteristics resulted from the cell growth during the fermentation process. The oxygen consumption of the microbial cells is considered the disturbance to the control system and is estimated by using the mathematical model and compensated by the FF control action. The substrate is FB controlled by repeated pulse-fed of carbon source. The schematic diagram for the control system is shown in Figure 5 [13].
The specific cell growth rate is modeled using double substrate Monod equation shown by Eq. (61). The equations for the specific glucose consumption rate, the specific oxygen consumption rate, OUR, and OTR are shown by Eqs. (62)-(65), which are the same as that of Section 4.2.1.
The mass balance equations for the repeated fed-batch culture are described by Eqs. (66-69).
where S F is the substrate concentration in the concentrated feeding solution. The substrate feeding solution is concentrated so that the volume change resulted from the substrate feeding can be neglected in Eqs. (66)-(68).
For FF control of DO, in order to compensate the DO disturbance resulted from the cell growth, OTR should be equal to OUR according to Eq. (68) if the dilution effect of the feeding is neglected so as to ensure C L unchanged (dC L /dt ¼ 0) and remained at the set-point. As the oxygen transfer driving force, ΔC ¼ (C L * À C L ), is relatively constant when C L is maintained at the set-point, k La should be controlled to meet Eq. (70) to compensate the time-varying OUR by the cell respiration according to Eqs. (65) and (68), and dC L /dt ¼ 0.
The value of k La is controlled by agitation speed (N) and aeration rate (G) shown by Eq. (71).
Between the two manipulated variables, N is more effective than G in controlling k La [14]. Therefore, 70% of the control effort is assigned to N and 30% is assigned to G by using Eqs. (72) and (73), respectively, which are drawn from Eqs. (70) and (71).
where the subscripts t and tÀ1 are the current and last time points, respectively. So, N and G control actions should be finished in several control rounds. Eqs. (72) and (73) are used in the FF control.
As model predictions may not be very accurate, FB control is used to eliminate the control error and ensure the control accuracy. In the case of FB control, the error between DO set-point and process variable is calculated by Eq. (74).
where C L.sp is DO set-point. The proportional and integration (PI) control strategy is used for FB control by using Eqs. (75) and (76) for N and G control, respectively. Similarly, 70% of the control action is assigned to N FB and 30% is assigned to G FB .
Then, the total DO control actions of N and G are shown by Eqs. (77) and (78) where N 0 and G 0 are the initial values of N and G, respectively.

Simulation results
In this system, DO is FF-FB controlled by agitation speed and aeration rate and the substrate concentration is FB controlled by repeated pulse-fed of certain amount of the concentrated feeding solution to make the substrate concentration reach 30 g/L when the substrate concentration is lower than the set-point of 5 g/L. In order to confirm the robustness of the control system under model prediction errors and noises, 5% randomized noises and 20% over estimate of the cell growth were added in the mathematical model predictions in FF control. Then, simulations were made with the above noises and prediction errors. The results indicated that even if the noises and relatively large model prediction errors existed, the control system still had good performance. The reason is that FB control finally compensated the inaccuracy of the FF control, as shown in Figure 6 [13].

Conclusion
Modeling and simulation are useful tools for understanding, analysis, and optimization of bioprocesses [14][15][16][17]. By using the modeling and computer simulation methods, the dynamics of cell growth and metabolism under different conditions and various fermenter operation modes can be evaluated and the information can be used for bioprocess optimization and bioreactor control.