Analysis of Hopf Points for a Zymomonas mobilis Continuous

Article pubs.acs.org/IECR

Analysis of Hopf Points for a Zymomonas mobilis Continuous Fermentation Process Producing Ethanol Hangzhou Wang, Nan Zhang, Tong Qiu, Jinsong Zhao, Xiaorong He, and Bingzhen Chen* Department of Chemical Engineering, Tsinghua University, 100084 Beijing, China ABSTRACT: Zymomonas mobilis is one of the most promising industrial microorganisms for industrial ethanol production. However, biomass, product and substrate oscillation behaviors have been observed under certain operating conditions in both experiments and numerical simulations of continuous fermentations involving this organism. This article analyzed oscillation behavior during fermentation based on Hopf bifurcation analysis, in which the Hopf singularity point and limit cycles were introduced. All Hopf points in the operating parameter values panel were simultaneously located and identified for given sets of operating conditions. Calculation of the distribution of Hopf points contributed to enhancing the stability of the fermentation process, maintaining high product quality, and providing insights into system dynamics.

■

minimizing ethanol inhibition.16 In addition, experimental and theoretical studies also illustrated that variables such as dilution rate and feed substrate concentration have an impact on static and dynamic behavior in continuous stirred tank fermentors.17,18 The existence of Hopf bifurcations were reported to cause oscillatory behavior in fermentation processes with Z. mobilis.2 Optimal operating points and regions, calculated from numerical analysis of oscillating biosystems, can be used to avoid the oscillatory behavior that unfavorably decreases the fermentation yield and declines product quality. Follow-up design and control strategies can be developed to maintain stability and guarantee high productivity and high product quality.6,18 With properly selected operating conditions, these Hopf points can be eliminated to avoid oscillations, directly leading to a more stable and productive fermentation process. The advantages of numerical analysis are listed as follows: 1. As the fermentation process is slow, its dynamic simulation through experiment may be inefficient and unable to locate points corresponding to certain dynamic behaviors, such as bifurcation and oscillation. 2. Because only a limited number of experiments can be performed, some dynamic behavior may be missed or neglected. Existing experimental work has shown the existence of oscillatory behavior; researchers analyzed these dynamic systems using single operation parameter analysis with all other state variable changes monitored,1,2,14,16,18,19 detecting frequent bifurcation, chaotic behavior, and even force oscillations.1,14 However, all previous work has suffered limitations because of the lack of ability to deal with two or three simultaneously changed operation parameters. In this article, the effect of operating condition values on product yield was analyzed, with

INTRODUCTION Bioethanol is an important renewable and sustainable alternative fuel source,1 both directly in the form of fuel ethanol and as a blend with gasoline. Agricultural feed reactants that are suitable for fermentation processes producing bioethanol contain sugars, starches, and cellulose materials. The yeast Saccharomyces cerevisiae2−5 and the bacterium Zymomonas mobilis (Z. mobilis)3−6 are currently the most important microorganisms for industrial production of ethanol. Z. mobilis produces less biomass than Saccharomyces cerevisiae, while more carbon is funneled to the ethanol fermentation. It was reported that the ethanol yield of Z. mobilis could be as high as 97% of the theoretical yield of ethanol to glucose,7 whereas only 90−93% can be achieved by Saccharomyces cerevisiae.6 Thus, Z. mobilis has been promoted as a more promising microorganism than yeast for the industrial production of ethanol.8 Continuous fermentation using Z. mobilis is prone to dynamic instability,2,3,6,9 in which sustained oscillations of biomass, substrate, and product concentrations are routinely observed, especially at low growth rate and high ethanol concentrations.10 Oscillations are mainly caused by product inhibition of enzymatic reactions, and experimental studies have revealed that high temperatures may enhance such inhibitory effects.11 Meanwhile, previous work has also showed that oscillatory states are achieved at low dilution rates.10 These oscillations increase the average residual glucose, which correspondingly decreases the ethanol yield and declines product quality.6 A considerable amount of work2,3,6,7,10−15 has demonstrated oscillatory behavior in fermentation processes both experimentally and numerically. With regard to bioethanol fermentation systems utilizing Z. mobilis, several experimental and theoretical studies focused on dynamic stability, which affects the industrial productivity of fuel ethanol.3 In alcoholic fermentation, it is assumed that ethanol inhibits the maximal specific growth rate indirectly by inhibiting the synthesis of the main intermediate, which in turn determines the overall growth rate.11 Therefore, the glucose conversion, the ethanol yield, and the production rate will be favored by © 2012 American Chemical Society

Received: Revised: Accepted: Published: 1645

May 17, 2012 December 12, 2012 December 18, 2012 December 18, 2012 dx.doi.org/10.1021/ie3013049 | Ind. Eng. Chem. Res. 2013, 52, 1645−1655

Industrial & Engineering Chemistry Research

Article

the oscillatory phenomenon is one of the most significant. In the oscillatory situation, state variables experience periodic changes. Theoretically, oscillation corresponds to Hopf points pertaining to singular points. A supercritical bifurcation occurs, with orbital limit cycles, when the reaction system is situated near this type of point. A pair of conjugate eigenvalues in mirror symmetry, with at least one nonzero real part in the other n−2 eigenvalues, can indicate this phenomenon. Periodic systematic oscillations will be detected near these Hopf singularity points, which clarifies the significance of locating Hopf points. Currently, most research on Hopf points emphasizes accurate search for and location of Hopf points and limit cycles.

rigorous Hopf points location in the operating region. Moreover, the synergism between three multivariate operating parameters was discussed. As stated, numerical analysis methods can enable complete and thorough analysis of dynamic systems. In this article, we first introduce the numerical methodology for locating Hopf point in steady state solutions as well as the limit cycles and then apply the proposed methodology to the Z. mobilis fermentation process.

■

INTRODUCTION TO HOPF POINTS AND NUMERICAL ANALYSIS METHODS In general, chemical processes are assumed to be described by the following equation: ⎧ dx = F (x , α ) ⎪ ⎨ dt ⎪ ⎩ x(0) = x 0

■

HOPF POINT DETECTION AND LOCATION The unique feature of a Hopf point lies in its special pair of conjugate pure image eigenvalues associated with the Jacobian matrix at the steady state point. By the Hopf bifurcation theorem, Hopf point detection can be achieved by calculating eigenvalues. A more common way is to design a detection function to identify Hopf points and the corresponding unique stable steady state. Consider the function described in eq 3

(1)

where x represents state variables in the dynamic system such as temperatures, pressures, components concentration, etc. (x ∈ Rn) and α is a parameter representing the operating condition value or design variable, such as flow rate, initial feedstock concentration, reactor volume, and so on. Equation 1 describes the dynamic behavior of a set of state variables under specified operating conditions, providing access to better understanding of system behaviors and optimization. The solution satisfying eq 2 is the equilibrium solution representing a set of practical state variable values under steady operating conditions.

F (x , α ) = 0

ψH(x , α) =

∏ (λi(x, α) + λj(x, α)) (3)

i>j

where λi are the eigenvalues of the Jacobian matrix at a specified steady state point. According to theoretical analysis, this function value will vanish when encountering a Hopf bifurcation point, associated with a pair of conjugate pure image eigenvalues, λ1,2 = ± iω0. Note that paired real eigenvalues λ1 = κ, λ2 = −κ can be a special case of ψH = 0, which corresponds to a neutral saddle situation, rather than the Hopf point we are primarily concerned with. In that sense, the Hopf point must be carefully identified from neutral saddle points when the zero point of detection function is found. At a generic Hopf bifurcation point, the function ψH is expected to be real, smooth, and have a regular zero value, with a nonsingular Jacobian of the proposed dynamic system assured simultaneously at this specific point. Newton’s method can be applied, and to avoid tedious explicit computation of all the eigenvalues of A, the current work used a method based on the bialternate product to compute ψH. The procedure of Hopf point detection based on the bialternate product is expressed in eq 4:

(2)

As the main operating points, a steady state solution set inevitably contains singular points that are worth looking into to enable smooth operation. A thorough study of singularity helps to screen and evade these distinct points, at which rich dynamic behavior exists and stability is practically difficult to maintain. Meanwhile, the accurate location of singular points, including the clear distribution of Hopf points, will no doubt contribute to safer decisions when a change in operating conditions must be made, because a number of specific singular points may possibly appear when switching across the operating region. Here, the matrix Fx(x,α), which is the Jacobian matrix of F(x,α), is introduced to assess stabilities of steady state operation points. If x = x* is the solution of F(x,α) = 0, then under the operating condition governed by α, the system will maintain in the steady state when all state variables take a full set values of x. In general, steady state points have different stability characteristics and the stability of a system under certain operating conditions and state variable values is determined by the eigenvalues of the Jacobian matrix. Specifically, the system is stable if all eigenvalues are located on the left of the complex panel, otherwise the system is unstable. In other words, when every real part of all the Jacobian matrix eigenvalues is negative, re(λi) < 0, i = 1,2, ...,n, then the system is stable when a minor disturbance is encountered. On the contrary, if there is only one eigenvalue that consists of a positive real part, the whole system would fall into an unstable situation under which even a small disturbance would cause extensive system changes. Stable steady state operating conditions are desirable in both operation and design, but some specific unstable steady situations also have potentially practical and theoretical applications, among which

ψH(x , α) = det(2Fx(x , α) ⊙ In)

(4)

Here, ⊙ represents the bialternate product of the matrix and formulated as apr aps (2A ⊙ In)(p , q),(r , s) =

δqr δqs

+

δpr δps aqr aqs

(5)

where ⎧1, i = j δij = ⎨ ⎩ 0, i ≠ j ⎪

⎪

In the matrix, rows are labeled by the multi-index (p, q) (p = 2, 3, ..., n; q = 1, 2, ..., p − 1), and columns are labeled by the multi-index (r, s) (r = 2, 3, ..., n; s = 1, 2, ..., r − 1). The result of the determinants is enumerated in eq 6: 1646

dx.doi.org/10.1021/ie3013049 | Ind. Eng. Chem. Res. 2013, 52, 1645−1655


Article

Figure 1. Schematic of Hopf analysis procedure.

(2A ⊙ In)(p , q),(r , s)

⎧−aps , ⎪ ⎪ apr , ⎪ ⎪ apr + aqs , =⎨ ⎪ aqs , ⎪ ⎪−aqr , ⎪ ⎩ 0,

■

r=q

LOCATION OF LIMIT CYCLES Bifurcation theory suggests that limit cycles are intrinsically of concomitance with Hopf singularity points. The shooting method, a widely accepted algorithm for calculating limit cycles, first converts the proposed question to a periodic boundaryvalue problem for the cycle continuation in eq 7 as follows:

r ≠ q, s = q r = p, s = q r = p, s = q p=s otherwise

⎧ dx T F (x , α ) = ⎪ 2π ⎨ dt ⎪ x(0) = η ⎩

(6)

Theoretically, the detection function equals zero when the Hopf bifurcation occurs, but in practice, we first identify the sign of the test function as positive or negative. If ψH(xi,αi)ψH(xi+1,αi+1) < 0, this reveals a Hopf singularity point standing between (xi,αi) and (xi+1,αi+1); thus, the rigorous Hopf point (xH,αH) location is subsequently found through a gradually reduced search space.

(7)

After parameter α is preset, we search for an initial value η and a period T satisfying eq 7 to obtain the periodic solution for the problem. Fortunately, software packages such as AUTO, 20−22 CONTENT,23 MATCONT,24,25 DsTool,26 PyDSTool,27 XPPAUT,28 and some previous research by the authors’ 1647



Article

Figure 2. Fermentation process.

group29 have implemented the method described and laid the foundations for the current work. The strategy for locating Hopf points in a dynamic system was summarized in Figure 1. Beginning with a model of the fermentation process, this procedure was successively conducted by calculating the steady state solution of the dynamic system as parameters were changed, then identifying the steady state solutions to determine whether or not there are Hopf points. If there is a Hopf point and the first Lyapunov coefficient is negative, orbital stable limit cycles in the system and corresponding oscillatory behavior are sure to exist in the dynamic system. Finally, parameters were varied to find all Hopf singularity points by this ergodic process, providing a valuable reference for operating point optimization. Fermentation Process Flowchart. There are several reported models10,11,13−15 describe the Z. mobilis fermentation process. Figure 2 shows a diagram of the continuous flow fermentor, while the data used in the numerical simulation is listed in Table 1. The kinetics developed by Jobses et al. for production of bioethanol by Z. mobilis was adopted.1 Product Inhibition. Ethanol production by Z. mobilis is subjected to end-product inhibition. The product alters the cell membrane composition and inhibits enzymatic reactions such as carrier-mediated transport processes and metabolic syn-

theses. In this work, to decrease inhibition by ethanol, a removal membrane was positioned at the bottom of the fermentor to remove ethanol produced during fermentation. Continuous Ethanol Removal. A coupled membrane separation system for ethanol removal was proposed to prevent ethanol inhibition and increase the productivity and efficiency of the fermentation process. The membrane has selectivity toward ethanol within the reaction environment, with a gas/ liquid sweep stream on the permeate side to remove product from the membrane surface. With the removal of ethanol, inhibition would be prevented and cell activity maintained at a high level, making it possible to use a higher glucose feed concentration, resulting in increased fermentor cell density. The increased cell activity and cell density contribute to higher volumetric productivity and a smaller fermentor volume. Additionally, the use of a concentrated glucose feed will decrease the amount of process water.19 Model Assumptions. The model adopted in this paper was based on the following simplified assumptions:30 1. The fermentor is considered to be a well-mixed system and is not controlled by diffusion. 2. The feed contains only nutrient (substrate). 3. The specific growth rate is defined by a Monod-type equation. 4. The maximum specific growth rate is proportional to the concentration of the key component (e). 5. k1, k2, k3, Ks, ms, mp, Ysx, and Ypx are constants. Model Development and Discussion. One of the most generally acknowledged models of fermentation processes is the maintenance model, for which substrate consumption equation is as follows:

Table 1. Parameters in the Fermentation Process Model param.

value

unit

Cs0 DMin AM Din Ks k1 k2 k3 ms mp Ysx Ypx Cp0 Ce0 Cx0 CpM0 ρ VM VF Pm

140 0.5 0.24 0.04 0.5 16 0.497 0.00383 2.16 1.1 0.0244498 0.0526315 0 0 0 0 789 0.0003 0.003 0.1283

kg/m3 h−1 m2 h−1 kg/m3 h−1 m3/kg·h m6/kg2·h kg/kg·h kg/kg·h kg/kg kg/kg kg/m3 kg/m3 kg/m3 kg/m3 kg/m3 m3 m3 m/h

⎛ μ ⎞ + ms⎟Cx rs = ⎜ ⎝ Ysx ⎠

where μ is the specific growth rate. The first term represents the growth rate and the second term represents the maintenance. Similarly, the product generation rate equation is defined as

⎛ μ ⎞ + mp⎟⎟Cx rp = ⎜⎜ ⎝ Ypx ⎠ The growth rate of biomass rx has its classical definition and can be expressed in the form:

rx = μCx The rate of the key compound re is formulated as 1648



Article

where the subscripts Din, Dout, DMin, and DMout represent their respective dilution rates, and Ysx and Ypx refer to the yield factors of biomass on substrate and product; k1, k2, k3 are empirical constants; V stands for volume; mS and mP are maintenance factors based on substrate requirements and product formation. The membrane side and the fermentation side are distinguished by subscripts M and F.

re = f (Cs)f (Cp)Ce

The formation rate of the key component is a function of the substrate concentration, inhibited by ethanol, while its activity is expressed in terms of concentration only.11 Subscript s represents the substrate, while subscript p refers to the product. f(Cs) is given by the Monod equation: f (Cs) =

■

Cs K s + Cs

HOPF POINTS IN THE SYSTEM AND THE LIMIT CYCLES The set of steady state points in the investigated process of fermentation is shown in Figure 3. In the case of Cs0 = 150.3

The function f(Cp) is defined in polynomial form: f (Cp) = k1 − k 2Cp + k 3Cp2

Thus, the equation for re can be substituted into the following expression: re =

k1 − k 2Cp + k 3Cp2 K s + Cs

CsCe

The base set of parameters in the current investigation is listed in Table 1.19 The dynamic model for the continuous fermentation process for substrate (s), microorganism (x), and the key compound (e) are given by the following set of ordinary differential equations: dCx = μCx + DinCx0 − Dout Cx dt ⎛ μ ⎞ dCs = −⎜ + ms⎟Cx + DinCs0 − Dout Cs dt ⎝ Ysx ⎠ Figure 3. Changes in product concentration with dilution rate.

k1 − k 2Cp + k 3Cp2 dCe = CsCe + DinCe0 − Dout Ce dt K s + Cs

kg/m3, a Hopf point was detected at Din = 0.052 h−1, which coincides exactly with reported literature values.2,16 At the Hopf point, we have a product concentration Cp = 58.050 kg/m3, which approaches the maximum point of product yield of Din = 0.046 h−1, Cp = 58.103 kg/m3. This fact highlights the importance of Hopf point investigation in the present dynamic system, because it is so close to the maximum conversion point that we seek. In Figure 3 the Hopf singularity point is marked as a red star, with the first Lyapunov coefficient l1 value of −0.00153. The negative sign means the point corresponds to a supercritical Hopf bifurcation and the limit cycle is orbitally stable. In Figure 3, the abscissa denotes the feedstock dilution rate, ranging from zero to a reasonable maximum. As the variable indicating the degree of feedstock volume diffusing into the reactor volume, the dilution rate has a relatively small value in practice. Therefore, only a range of the dilution rate worthy of consideration is discussed. As mentioned, the Hopf points are theoretically very close to a risky but high conversion point, and associated limit cycles are depicted in Figure 4. Each limit cycle is generated by a very small disturbance from the initial operating parameter value at the Hopf point. According to bifurcation theory, experimentally observed or mathematically simulated oscillations in a dynamic system must start and end at certain critical points in the system. It has been proven that static and dynamic bifurcation behaviors, over a wide range of operating parameters, abound in the Z. mobilis fermentation process. In the current work, the fermentor dilution rate, the membrane dilution rate, the initial substrate

The product (ethanol) on the fermentor side (p), the product (ethanol) on the membrane permeate side (pM), and the mass balance for the product (ethanol) over the fermentor are given by ⎛ μ ⎞ = ⎜⎜ + mp⎟⎟Cx + DinCp0 − Dout Cp dt ⎝ Ypx ⎠ A P − M m (Cp − CpM) VF

dCp

The mass balance for the product (ethanol) over the membrane is given by dCpM dt

=

AM Pm (Cp − CpM) + DMinCpM0 − DMout CpM VF

where the AM is the permeation area and Pm is the permeability of the membrane for ethanol. The exit fermentor dilution rate can be obtained as follows: Dout = Din −

AM Pm (Cp − CpM) VFρ

The corresponding exit membrane dilution rate becomes DMout = DMin +

AM Pm (Cp − CpM) VMρ 1649



Article

potential oscillatory operation points can also be measured to better estimate the operating risk as a probability. Case 1: Cs0 and Din. In this case, the dilution rate of the feedstock and the initial substrate concentration were examined because of their decisive roles as operating variables and their ease of manipulation. The steady state solutions were calculated; the Hopf points were marked in Figure 5, with

Figure 4. Homology limit cycles adjoining the Hopf point.

(glucose) concentration of feedstock concentration, and permeation area were chosen as bifurcation parameters, not only for their decisive role as operating variables in this system but also because of their ease of manipulation and control during practical fermentations. To avoid the oscillatory behaviors caused by the existence of the Hopf bifurcation, which as stated decrease the fermentation yield, system stability, and controllability, decline product quality, we put forward a rational method to position the entire set of Hopf points in the operation parameter panel. In this proposed method, all Hopf singularity points were located and assembled in the operation parameter panel as a “risky curve”. Thus, operating parameters can be located away from the undesirable singular points, unwitting steps across these points can be avoided, and appropriate optimized values of the parameters adopted.

Figure 5. Steady state points and Hopf points.

the first Lyapunov coefficient labeled to imply that these Hopf points are supercritical and will stimulate nearby periodic limit cycles. To observe the membrane-side dilution influence on the system by direct-viewing, results were combined within one bifurcation map, as shown in Figure 5, under variable membrane side dilution rates. In Figure 5, it can be seen that the membrane side dilution rate DMin changed from 0.1 to 0.9 h−1 at a step size of 0.2 h−1. Under each given membrane side dilution rate, production concentration is plotted against the feedstock dilution in the steady state solution. Note that the clusters of production concentration curves have generally delayed peaks with equal changes in membrane side dilution. As may be expected, the corresponding Hopf points were detected within the region of low feedstock flow rate, where the negative first Lyapunov coefficients indicate potential oscillation near the production concentration peak. Another fact worth noticing is that the product conversion rate peaks in the system decreased with increasing membrane side dilution rate, which can be explained by the conversion inhibition effect brought by increasing membrane side dilution rate. Fortunately, with the decrease in concentration of product, the Hopf point becomes increasingly separated from the conversion peak, characterized by an increasing first Lyapunov coefficient at the same time. Details of this case study and its five Hopf points are listed in Table 2. The same conclusion can be drawn after analyzing homologous state and operating variables. Figure 5 only describes the influence of a single parameter (or two parameters to a limited extent) on the Hopf points. However, it is more common for operating parameters to work together and change simultaneously. Thus, we examined

■

PARAMETER INFLUENCES ON HOPF POINTS Hopf points in Z. mobilis fermentation result in orbital stable limit cycles and oscillation behaviors in the process, and it follows from this, as mentioned, that the product conversion rate, system stability, controllability, and product quality are all adversely affected. Here, four important operating parameters, including dilution rate of feedstock, membrane side feedstock dilution rate, the initial substrate concentration, and the permeation area, were investigated. Each case study was undertaken in the following way. First, the changes in product concentration were observed and calculated under various operating conditions and the results shown in a single figure to confirm the existence of the Hopf points along with the steady state operation points in the system. These two selected parameters were then changed over a range wide enough to find all Hopf points in the operating region formed orthogonally by these two parameters. As can be seen, every point in the figure was obtained under fixed corresponding operating parameters. Among all steady state solutions of the system, a Hopf point intrinsically exists with joint orbital stable limit cycles in practice, arising from the oscillation. Furthermore, using the methodology proposed in this paper, not only can these undesirable points be avoided but the distances between current operating conditions and 1650



Article

Table 2. Hopf Points in the System Cs0 = 140 kg/m3, AM = 0.24 m2 DMin (h−1)

Din (h−1)

Cp (kg/m3)

CpM (kg/m3)

first Lyapunov coefficient l1(10−4)

0.1 0.3 0.5 0.7 0.9

0.096 0.395 0.627 0.832 1.027

57.228 53.154 50.790 48.995 47.469

57.168 52.988 50.527 48.641 47.030

−8.968 −2.583 −1.825 −1.344 −1.073

results in Figure 7 indicate that the membrane permeation area has only a slight influence on product concentration as a function of feedstock dilution rate.

another important parameter, the initial substrate concentration. Specifically, while we fixed AM = 0.24 m2 and varied DMin from 0.1 to 0.9 with a step-size of 0.2 h−1, we searched Din and Cs0 to find the existence of Hopf points and then recorded these latter two parameter values to plot the data in one figure. Accordingly, all parameter values creating Hopf singularity points over the range of fixed DMin and AM values are shown in Figure 6.


As can be seen in Figure 7, although the permeation areas differed markedly, their impact on product concentration as a function of feedstock dilution was small and the steady state solution curves almost coincided. In a similar manner, the calculated Hopf points were located very close to one another and were almost constant. Thus, it can be concluded that the permeation area was not a sensitive parameter for this system. From Table 3, it can be seen that AM contributes little to Cp and CpM, as shown in Figure 7. Table 3. Hopf Points in the System with Cs0 = 140 kg/m3, DMin = 0.5 h−1 Figure 6. All Hopf points in the parameter panel.

In Figure 6, the Hopf points can be seen to change in response to the variations in feedstock dilution rate and initial substrate concentration. When the feedstock dilution rate or the initial substrate concentration was low, the system displayed a higher probability of generating Hopf points. Notably, even if one parameter was altered only over a certain narrow range, while the other was varied more aggressively, Hopf points were located, which reveals that there is no universal safe operating region for every parameter. The position of the Hopf points, as expected, changed with dilution rate on the membrane side. It is demonstrated in Figure 6 that a greater membrane side dilution rate causes a heavier density of Hopf points within the operating range as the feedstock dilution rate increases. Meanwhile, locations across the Hopf points are shown in this figure by a dashed line, which marks the Hopf points under a constant substrate concentration of 140 kg/m3. After examining D Min , the final parameter A M was investigated because the diversity of membranes used in fermentation processes may cause different permeation areas of the membrane. Here, we selected AM from 0.1 to 0.9 at intervals of 0.2 m2, and fixed DMin = 0.5 h−1, to examine the influence of membrane area on product concentration. The

AM (m2)

Din (h−1)

Cp (kg/m3)

CpM (kg/m3)


0.1 0.3 0.5 0.7 0.9

0.650 0.624 0.620 0.618 0.617

50.578 50.817 50.861 50.879 50.889

49.954 50.607 50.734 50.788 50.818

−1.564 −1.770 −1.825 −1.847 −1.827

Furthermore, at a fixed value of permeation area, Hopf points were fixed within the operating panel with respect to initial substrate concentration and feedstock dilution rate, as shown in both Figures 7 and 8. Thus, in the analysis that follows, the permeation area was not considered because it is not a sensitive variable in this system. In contrast, the dilution rate on both the feedstock and the membrane sides, together with the initial substrate concentration in the feed were discussed in more detail. Case 2: DMin and Cs0. In case 2, we focused on membrane side dilution and the initial substrate concentration in the feedstock. Simulation was carried out under five given feedstock dilution rates and one membrane-side dilution rate. The product concentration changed with the initial feedstock substrate concentration, as shown in Figure 9. 1651



Article

Table 4. Hopf Points in the System AM = 0.24 m2, DMin = 0.5 h−1 Din (h−1)

Cs0 (kg/m3)

Cp (kg/m3)

CpM (kg/m3)


0.02 0.05 0.1 0.3 0.5

418.671 250.692 195.376 155.680 144.557

58.833 58.201 57.267 54.330 52.050

58.526 57.897 56.968 54.047 51.780

−3.048 −1.375 −1.070 −1.222 −1.574

Figure 8. All Hopf points in the parameter panel.


From Figure 10, it can be seen that, once again, the entire system has a turning point on each curve. The position of the Hopf points, although always below the turning point, changes markedly as the initial substrate concentration increases, whereas the rates decreased above the turning points. The dilution rate of feedstock influenced the curve shape to a large extent, with larger dilution values giving lower initial substrate concentrations at which Hopf points were encountered. Case 3: Din and DMin. In this case, the dilution rate on both sides of the membrane were considered. Calculation results are depicted in Figure 11, showing a complete three-parameter bifurcation diagram. When the initial substrate concentration and the feedstock dilution rate were both fixed, the product concentration decreased progressively with the membrane side dilution rate because the feedstock dilution rate had only a low value (0.02 h−1), which gave a sufficiently long resident time for the product concentration to achieve its maximum. With an elevated membrane-side dilution rate, the ethanol was gradually removed through entrainment, resulting in a drop of the product concentration. This effect would be intensified by further increases in the membrane-side dilution rate. Meanwhile, when the initial substrate concentration in the feedstock increased, the product could maintain maximum growth rate for a longer time, even at the same membrane-side dilution rate. All Hopf points are marked in Figure 11, and actual values are listed in Table 5. Contrary to to the preceding cases, Table 5 clearly indicates that Hopf points are generated at the maximum growth rate, adjacent to the value of product concentration, as for the membrane-side product concentration. The Hopf points in the dilution rate and membrane-side dilution rate panel are


Product concentration increased gradually with the increase of initial substrate concentration of feedstock. And the turning points seem to serve as an upper limit of concentration no matter how the substrate concentration increases. This baffle effect can be interpreted since at the begin in the fermentation process, the Z. mobilis grow with the increase of initial concentration of feedstock, but when the growth reaches the maximum, the produce capability has no potential to constantly increase. However, the dilution rate of feedstock also has an impact on product concentration. It is worth noting that the greater feedstock dilution rate is, the earlier turning points will appear, with low concentration at high dilution rate. Hopf points, shown in Figure 9, tighten therewith the maximum point. Corresponding values are listed in Table 4. Similarly, in Table 4, we find all the first Lyapunov coefficients are negative, indicating the existence of stable limit cycles. The location of Hopf singularity points in the parameter panel changed significantly with the initial substrate concentration and the membrane side dilution rate, as shown in Figure 10. 1652



Article

feedstock and membrane side, together with the initial substrate concentration, are the most sensitive parameters contributing to the product yield and product quality, and they significantly influence the location of Hopf points.

■

RESULTS AND DISCUSSION In our research work, four important operating parameters, including dilution rate of feedstock, membrane side feedstock dilution rate, the initial substrate concentration, and the permeation area, were investigated, because of their decisive roles as operating variables and their ease of manipulation. As operating conditions changed, the steady state solution curve and the location of Hopf point appear in very distinguishing manners. Specifically, we found that 1. Hopf points, in the present fermentation dynamic system, locate very close to the maximum conversion point, and this fact highlights the importance of Hopf point investigation. 2. Increasing the membrane-side dilution rate will result in delayed and decreased the product conversion rate peak value, while the undesired Hopf point becomes increasingly separated from the desired conversion peak. 3. Hopf points have a high probability to appear in areas of either low feedstock dilution rate or low initial substrate concentration; thus, there is no universal safe operating region for a single parameter. 4. The membrane permeation area is not a sensitive parameter and only has a slight influence on product concentration. Almost coincident steady state solution curves and Hopf points were detected. 5. Baffle effect of product concentration and the maximum turning point were observed, while the position of Hopf points changes markedly as initial substrate concentration increases. 6. The three-parameter bifurcation diagram showed that when the initial substrate concentration is fixed and the feedstock dilution rate is low, the product concentration decreased progressively with the membrane side dilution rate, but when the initial substrate concentration in the feedstock increased, the product could maintain a maximum growth rate for a longer time. 7. At low initial concentration, the Hopf points appear at almost the same level of feedstock dilution rate, whereas at higher initial concentrations, the Hopf points arise in a region of relatively low dilution rate but high membrane dilution. 8. Product conversion rate, system stability, controllability, and product quality all depend heavily on these four operating variables. The synergism between multivariate operating parameters has been first discussed. The emphasis of this article is on locating Hopf singularity points in the process. As a result, with four major manipulate parameters and their combinations, also within a large operation condition parameters scopes, we have studied the Hopf singularity points thoroughly in our research work. The results showed the universal existence of Hopf singularity points in this process within a large operational space. Although it seems there were many singularity points in the panel constructed by operational variable, suitable operation points for designing process were even more in the panel. Because only the points located exactly in the critical curve were not recommended, many more points away from the curve in the


Table 5. Hopf Points in the System Din = 0.02 h−1, AM = 0.24 m2 Cs0 (kg/m3)

DMin (h−1)

Cp (kg/m3)

CpM (kg/m3)


150 200 250 300 350

0.049 0.133 0.218 0.302 0.386

58.701 58.731 58.757 58.781 58.803

58.671 58.648 58.622 58.594 58.565

−5.710 −3.884 −2.598 −1.651 −0.949

illustrated in Figure 12. Hopf points exist at each given initial substrate concentration with various dilutions.


It can be observed from Figure 12 that at low initial concentration, the Hopf points appear at almost the same level of feedstock dilution rate, whereas at higher initial concentrations, the Hopf points arise in a region of relatively low dilution rate but high membrane dilution. From all of these analyses, we can now conclude that, among all operating parameters, the dilution rates on both the 1653



Article

Notes

panel are feasible to the process design. Also, in practice, the operation condition parameter space was limited into a smaller operational space, which would meets complex constraints; the operation variables and their combinations also will be reduced, so there would be less Hopf singularity points, even no Hopf singularity points would be encountered in some situations. In general situations at process design phase, it is a better method to determine the operation ranges first, after that, to use the method in this article to detect whether there are Hopf singularity points or not. Most of the time, there would be seldom Hopf singularity points within the small feasible operational range that meets all the complicated practical constraints. Even in the worst situation, there was a Hopf singularity curve in the operational space; these points are feasible to process design except that these operation points just locate in the critical curve; for insurance, these points very close to the critical curve also can be excluded in process design.

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS The authors gratefully acknowledge support from China Postdoctoral Science Foundation (No. 20100480282), China Postdoctoral Science special Foundation (No. 2012T50099), and the National Basic Research Programme (No. 2012CB720500).

■

■

CONCLUSION In this article, a typical fermentation process using Z. mobilis to produce ethanol was investigated. Although it is a more promising microorganism than yeast for the industrial production of ethanol, continuous fermentation of Z. mobilis undergoes oscillatory behavior in which biomass, product, and substrate cycle under certain fermentation conditions. Therefore, it is important to search for a reasonable approach that helps to avoid the oscillatory behaviors caused by Hopf bifurcations, which decrease fermentation yield, system stability and controllability, and decline product quality. This article rigorously analyzed the Hopf point variation under carefully designed operating conditions and discussed its influence on product yield and quality. The article primarily focused on the operating conditions and product yield rather than complicated theoretical bifurcations and chaotic behavior. A numerical analysis method for locating the Hopf singularity point and limit cycles was introduced and applied to the Z. mobilis fermentation process for producing ethanol. In our research, the entire set of Hopf points have been solved in the operating parameter panel. Furthermore, to avoid tedious explicit computation of all the eigenvalues of a matrix, the current work used a method based on the bialternate product to compute the test function to determine Hopf point. It is the proposed method makes it possible for us to first locate and assemble these Hopf points as a “risk curve”, with reduced computational work and better accuracy. The calculated distribution of Hopf points can enable all operating parameters to be located away from undesirable singular points and avoidance of unwitting steps across these points, with appropriate and optimized values. Furthermore, the proposed method can be applied to other fermentation processes and improved models. Design and control strategies can be developed at a very early stage, due to the a priori identification of potential Hopf points, and this will contribute to the achievement of high productivity, product quality, and stability in dynamic processes.

■

NOMENCLATURE Ci Concentration of component i (kg/m3) ri Production rate of component i (kg/m3 h) Cs0 Inlet substrate concentration (kg/m3) DMin Inlet membrane dilution rate, inlet flow rate/membrane volume (h−1) DMout Output membrane dilution rate, output flow rate/ membrane volume (h−1) Din Inlet fermentor dilution rate, inlet flow rate/for monitor volume (h−1) Dout Output fermentor dilution rate, output flow rate/ fermentor volume (h−1) AM Area of membrane (m2) XS Substrate conversion KS Monod constant (kg/m3) k1 Empirical constant (h−1) k2 Empirical constant (m3/kg h) k3 Empirical constant (m6/kg2 h) ms Maintenance factor based on substrate (kg/kg·h) mp Maintenance factor based on product (kg/kg·h) Ysx Yield factor based on substrate (kg/kg) Ypx Yield factor based on product (kg/kg) VM Membrane volume (m3) VF Fermentor volume (m3) Pm Membrane permeability (m/h)

Greek Symbols

ρ Ethanol density (kg/m3) μ Specific growth rate (h−1) Subscripts i

■

s Substrate (glucose) inside the fermentor s0 Influent glucose to the fermentor p Product (ethanol) inside the fermentor p0 Influent ethanol to the fermentor e Internal key component inside the fermentor e0 Influent internal key component to the fermentor x Biomass (microorganisms) inside the fermentor x0 Influent biomass to the fermentor pM Product (ethanol) inside the membrane pM0 Influent ethanol to the membrane

REFERENCES

(1) Abashar, M. E. E.; Elnashaie, S. S. E. H. Multistablity, bistability, and bubbles phenomena in a periodically forced ethanol fermentor. Chem. Eng. Sci. 2011, 66 (23), 6146−6158. (2) Sridhar, L. N. Elimination of oscillations in fermentation processes. AIChE J. 2011, 57 (9), 2397−2405. (3) Paz Astudillo, I. C.; Cardona Alzate, C. A. Importance of stability study of continuous systems for ethanol production. J. Biotechnol. 2011, 151 (1), 43−55. (4) Shen, Y.; Zhao, X. Q.; Ge, X. M.; Bai, F. W. Metabolic flux and cell cycle analysis indicating new mechanism underlying process oscillation in continuous ethanol fermentation with Saccharomyces

AUTHOR INFORMATION

Corresponding Author

*Tel.: +86 10 62784572. Fax: +86 10 62770304. E-mail: [email protected]. 1654



Article

cerevisiae under VHG conditions. Biotechnol. Adv. 2009, 27 (6), 1118− 1123. (5) Bai, F. W.; Ge, X. M.; Anderson, W. A.; Moo-Young, M. Parameter oscillation attenuation and mechanism exploration for continuous VHG ethanol fermentation. Biotechnol. Bioeng. 2009, 102 (1), 113−121. (6) Bai, F. W.; Anderson, W. A.; Moo-Young, M. Ethanol fermentation technologies from sugar and starch feedstocks. Biotechnol. Adv. 2008, 26 (1), 89−105. (7) Sprenger, G. A. Carbohydrate metabolism in Zymomonas mobilis: A catabolic highway with some scenic routes. FEMS Microbiol. Lett. 1996, 145 (3), 301−307. (8) Li, C. C. Mathematical models of ethanol inhibition effects during alcohol fermentation. Nonlinear Anal.: Theory, Methods, Appl. 2009, 71 (12), 1608−1619. (9) Diehl, F. C.; Trierweiler, J. O. Control strategy for a Zymomonas mobilis bioreactor used in ethanol production. In Computer Aided Chemical Engineering; de Brito Alves, R. M.; Oller do Nascimento, C. A.; Biscaia, E. C., Eds.; Elsevier: Amsterdam, 2009; Vol. 27, pp 1605− 1610. (10) Ghommidh, C.; Vaija, J.; Bolarinwa, S.; Navarro, J. M. Oscillatory behavior of Zymomonas in continuous cultures: A simple stochastic model. Biotechnol. Lett. 1989, 11 (9), 659−664. (11) Jöbses, I.; Egberts, G.; Luyben, K.; Roels, J. Fermentation kinetics of Zymomonas mobilis at high ethanol concentrations: Oscillations in continuous cultures. Biotechnol. Bioeng. 1986, 28 (6), 868−877. (12) Bai, F. W.; Chen, L. J.; Anderson, W. A.; Moo-Young, M. Parameter oscillations in a very high gravity medium continuous ethanol fermentation and their attenuation on a multistage packed column bioreactor system. Biotechnol. Bioeng. 2004, 88 (5), 558−566. (13) McLellan, P. J.; Daugulis, A. J.; Li, J. H. The incidence of oscillatory behavior in the continuous fermentation of Zymomonas mobilis. Biotechnol. Prog. 1999, 15 (4), 667−680. (14) Daugulis, A. J.; McLellan, P. J.; Li, J. H. Experimental investigation and modeling of oscillatory behavior in the continuous culture of Zymomonas mobilis. Biotechnol. Bioeng. 1997, 56 (1), 99− 105. (15) Jarzebski, A. B. Modeling of oscillatory behavior in continuous ethanol fermentation. Biotechnol. Lett. 1992, 14 (2), 137−142. (16) Garhyan, P.; Elnashaie, S. Static/dynamic bifurcation and chaotic behavior of an ethanol fermentor. Ind. Eng. Chem. Res. 2004, 43 (5), 1260−1273. (17) Garhyan, P.; Elnashaie, S. S. E. H.; Al-Haddad, S. M.; Ibrahim, G.; Elshishini, S. S. Exploration and exploitation of bifurcation/chaotic behavior of a continuous fermentor for the production of ethanol. Chem. Eng. Sci. 2003, 58 (8), 1479−1496. (18) Garhyan, P.; Elnashaie, S. S. E. H. Experimental investigation and confirmation of static/dynamic bifurcation behavior in a continuous ethanol fermentor. practical relevance of bifurcation and the contribution of Harmon Ray. Ind. Eng. Chem. Res. 2004, 44 (8), 2525−2531. (19) Mahecha-Botero, A.; Garhyan, P.; Elnashaie, S. Non-linear characteristics of a membrane fermentor for ethanol production and their implications. Nonlinear Anal.: Real World Appl. 2006, 7 (3), 432− 457. (20) Doedel, E. J. AUTO: A program for the automatic bifurcation analysis of autonomous systems. Congr. Numer. 1981, 30, 265−284. (21) Doedel, E.; Kernevez, J. P. AUTO: Software for continuation and bifurcation problems in ordinary differential equations, including the AUTO 86 user manual. In Applied Mathematics; California Institute of Technology: Pasadena, CA, 1986. (22) Doedel, E.; Champneys, A.; Fairgrieve, T.; Kuznetsov, Y.; Sandstede, B.; Wang, X. AUTO 97: Continuation and Bifurcation Software for Ordinary Differential Equations, User’s Manual; Center for Research on Parallel Computing, California Institute of Technology: Pasadena, CA, 1997. (23) Kuznetsov, Y.; Levitin, V. CONTENT: Dynamical System Software; CAIN Europe: The Netherlands, 1997; available online:

http://www.computeralgebra.nl/systemsoverview/special/diffeqns/ content/gcbody.html. (24) Dhooge, A.; Govaerts, W.; Kuznetsov, Y. A. MATCONT: A MATLAB package for numerical bifurcation analysis of ODEs. ACM Transactions on Mathematical Software (TOMS) 2003, 29 (2), 141− 164. (25) Dhooge, A.; Govaerts, W.; Kuznetsov, Y. A. MATCONT: A Matlab package for numerical bifurcation analysis of ODEs. ACM SIGSAM Bull. 2004, 38 (1), 21−22. (26) Back, A.; Guckenheimer, J.; Myers, M.; Wicklin, F.; Worfolk, P. DsTool: Computer assisted exploration of dynamical systems. Notices Amer. Math. Soc 1992, 39 (4), 303−309. (27) Clewley, R. H.; Sherwood, W. E.; LaMar, M. D.; Guckenheimer, J. M. PyDSTool, A Software Environment for Dynamical Systems Modeling; 2007. http://pydstool.sourceforge.net (accessed 2012-5-9). (28) Ermentrout, B. Simulating, Analyzing, and Animating Dynamical Systems: A Guide to XPPAUT for Researchers and Students; Society for Industrial Mathematics: Philadelphia, PA, 2002; Vol. 14. (29) Wang, H. Z.; Chen, B. Z.; He, X. R.; Zhao, J. S.; Qiu, T. Numerical analysis tool for obtaining steady-state solutions and analyzing their stability characteristics for nonlinear dynamic systems. J. Chem. Eng. Jpn. 2010, 43 (4), 394−400. (30) Abashar, M.; Elnashaie, S. Dynamic and chaotic behavior of periodically forced fermentors for bioethanol production. Chem. Eng. Sci. 2010, 65 (16), 4894−4905.

1655


Analysis of Hopf Points for a Zymomonas mobilis Continuous

Recommend Documents