4314
Ind. Eng. Chem. Res. 1998, 37, 4314-4321
PROCESS DESIGN AND CONTROL Optimization of a Biochemical Fed-Batch ReactorsTransition from a Nonsingular to a Singular Problem Alok Jayant and S. Pushpavanam* Department of Chemical Engineering, Indian Institute of Technology, Kanpur 208016, India
In this work we analyze the operation of a biochemical reactor in the nonrepeated fed-batch mode. The reaction is assumed to follow Haldane kinetics (i.e., it is characterized by substrate inhibition). The feed rate of the substrate is chosen as the control variable. We optimize the system performance using Pontryagins’ maximum principle to obtain the flow rate profile F(t) for a nonlinear objective function. This is done without imposing any constraints on F(t) for the following cases: (1) fixed final volume, fixed final time; (2) free final volume, fixed final time; (3) fixed final volume, free final time; and (4) free final time, free final volume. The objective here is to maximize the biomass production such that the deviation from an average flow rate is minimal. This renders the problem nonsingular. The results of the optimal profile are compared with a semibatch reactor when a constant flow rate is maintained and also with the performance of a batch reactor. The effect of initial conditions on the different shapes of the optimal profiles is also discussed. Introduction Biochemical fermentation reactions are frequently characterized by substrate inhibition and/or product inhibition. The performance of these reactors is optimal when operated in the fed-batch mode, where the substrate is added to the reactor during the course of the reaction. The biomass is usually present initially in the reactor. The control variable is the substrate feed rate. The optimal performance of fermentation reactors in the fed-batch mode has been investigated with the objective of maximizing either (i) the biomass production or (ii) the product formed (Cazzador (1988), Cuthrell and Biegler (1989)). Two classical methods that have been used to study the problem are (i) Pontryagins’ continuous maximum principle and (ii) nonlinear optimization techniques as sequential quadratic programming (SQP) (Lee and Parulekar (1993), Levin (1992)). The model system most extensively investigated is characterized by a single reaction. The reaction rate varies linearly with biomass and is assumed to have a Haldane dependency on the substrate. The manipulated variable (i.e., the substrate feed rate) occurs linearly in the governing ordinary differential equations, which renders the problem singular. Hence, the optimal feed policy is either bang-bang or follows a singular arc, and the profile of manipulated variable can be discontinuous. The determination of the optimal policy requires solving different sets of equations along different regions as determined by the * Author to whom correspondence should be addressed. Current address: Department of Chemical Engineering, IIT, Madras 600036, India. E-mail:
[email protected]. Fax: (044) 2350509.
necessary and sufficient conditions of the maximum principle. Also, the numerical scheme is very sensitive to the initial guesses. In most problems, one almost needs to know the values of the variables being iterated upon before one can find them while using this procedure. Nonlinear programming techniques like SQP are based on discrete approximations of a continuous variable. Physical insight is used to formulate the problem and to render it amenable for solution using this technique. Here, the numerical method is computationally intensive and the problem with sensitivity to the initial guesses is eliminated. Weigand (1981) considered the problem of maximizing cell production in a repeated fed-batch mode of operation. He obtained an analytical expression relating the total reaction time with biomass concentration under the optimal conditions. In his analysis, however, he did not assume an upper bound on the permissible flow rate, which has led to his results being of only theoretical significance. Bonte et al. (1996) and Modak et al. (1986) determined the optimal feeding profiles for different kinetic models. They developed a computational scheme to obtain the feeding policy, that maximized (a) penicillin production and (b) cell mass production. The optimal control profile problem considered in these works is linear in the control variable, the substrate feed rate; hence, it is singular in nature. The optimal policy consists of different portions (i) when the feed rate is held at its minimum; (ii) when the feed rate is held at its maximum; (iii) when the feed rate follows a singular arc (Cazzador (1988)). Cazzador qualitatively analyzed the effect of initial conditions and substrate feed con-
10.1021/ie970557w CCC: $15.00 © 1998 American Chemical Society Published on Web 10/07/1998
Ind. Eng. Chem. Res., Vol. 37, No. 11, 1998 4315
centrations on the different profiles. His method of analysis can be applied only for one-dimensional problems, where the governing equations describe the evolution of a single scalar dependent variable. The numerical computation of the optimal profile for most problems, however, needs physical insight into the problem. This insight is necessary to determine the exact sequence of the different portions that constitute the optimal profile (Modak, 1986). Such an approach may not be applicable to systems where the interactions are very complex. In a nonbiochemical engineering context, optimal feed rate policies were obtained using pulse feeds at discrete intervals of time by Levin (1992). Sequential quadratic programming was used by Biegler (1984), Morrison and Sargent (1984), and Vasiliadis (1984). Shukla et al. (1998) used a sequential quadratic programming technique to obtain the optimum feed rate profile for a repeated fed-batch mode of operation. They divided the entire interval of operation into different segments. The feed rate in each interval was assumed to be a constant. The optimal feed rate policy thus determined does not use any information from the physics. This approach provides a general mathematical basis for obtaining the optimal solution. It is useful when the interactions between the different variables is complex and when it is not possible to infer the sequence of control policies that will give rise to optimal behavior. The disadvantage of this technique lies in the fact that it is computationally intensive. In this work we propose an alternative method to investigate the singular optimization problem discussed. The objective function is modified to include a term that is representative of the cost involved in preferring the optimal policy to the average constant flow rate; this is quadratic in the control variable. The problem is now nonsingular and can be analyzed elegantly using the maximum principle. This idea is based on concepts from quadratic control performance criteria, where a similar term is introduced in the design of controllers (Isermann, 1992). The classical singular problem investigated can be obtained as a limiting case of a nonsingular problem by reducing the weightage assigned to the quadratic term. This technique can be used to provide insight into the different regimes of the control variable and the sequence in which they occur. This information can then be used to solve the classical problem using the computational technique presented in Modak (1986). Thus, the motivation of this paper is to demonstrate an approach that will give us qualitative insight into the optimal feeding policy. This insight is then extrapolated to infer the manipulated variable profile in the limit of the singular problem. In this study we analyze the nonrepeated fed-batch mode of operation. We analyze four cases of reactor operation: (I) fixed final volume, fixed final time; (II) fixed final time, free final volume; (III) free final time, fixed final volume; and (IV) free final time, free final volume. We also discuss the effect of initial conditions on the optimal control profiles. Throughout this study we assume that there are no bounds on the control variable. Mathematical Model The evolution equations, which determine the reactor performance, are
(VX˙ ) ) µ(S)XV (VS˙ ) )
(1a)
-µ(S)XV + FSF Y
(1b)
(V˙ ) ) F
(1c)
Here, µ(S) represents the reaction rate dependency on the substrate concentration S, Y is the constant yield, SF is the substrate concentration in the feed, and F is the flow rate of the substrate. Here we assume
µ(S) )
µmaxS 2
RS + βS + γ
From eq 1a-c it follows that
(
)
(X˙ V) + (S˙ V) - (V˙ SF) ) 0 Y
(2a)
We now obtain
S ) SF +
(
)
V0 X0 X + S0 - SF - ) g(X,V) (2b) V Y Y
Here, the subscript ‘0’ is used to indicate the initial conditions in the reactor. This relationship is valid for all time and can be used to relate S(t) to X(t), which can be used to eliminate ‘S’ in the eqs 1a-c. Substituting in eqs 1a-c, we can reduce the three-dimensional system to the following two-dimensional system
X˙ ) Xµ (g(X,V)) -
FX V
V˙ ) F
(3a) (3b)
Our objective is to determine the feed rate profile F(t), that maximizes the biomass productivity at the terminal state. This result gives rise to the objective function
Minimize J ) - XfVf F(t)
(4a)
This problem however is singular because the control variable ‘F’ occurs linearly in eqs 3a and b. In this work we modify the objective function as follows
∫0t (F - F)2 dt
Minimize J ) - XfVf + F F(t)
f
(4b)
The subscript ‘f’ is used to denote the final state of the system. The introduction of the quadratic term in eq 4b renders the problem nonsingular. It can be interpreted as being representative of the cost involved in maintaining the optimal profile F(t) instead of the average profile F h . It is hence proportional to the deviation (F - F h )2. The scalar F is a weightage factor that can be used to change the importance of the two terms. The optimal control problem can be solved elegantly by using Pontryagins’ maximum principle. We define the Hamiltonian
[
H ) λX Xµ(g(X,V)) -
FX + λVF + F(F - F)2 (5) V
]
Here, λX and λV are the adjoint variables associated with eq 3a and b. Their evolution is governed by
4316 Ind. Eng. Chem. Res., Vol. 37, No. 11, 1998
[
∂µ(g(X,V)) F ∂X V
λ˙ X ) -λX µ(g(X,V)) + X
[
]
]
∂µ(g(X,V)) FX λ˙ V ) -λX X + 2 ∂V V
Un+1 ) Uni - J-1 i i Fi
(6a) (6b)
to generate the subsequent iterates. Here the values of Ui, Ji for the respective cases are
UI ) UII ) [λX0, λV0]t
The optimal flow rate profile is obtained by setting
H˙ F ) 0
∂λX ∂λX0 JI ) ∂V ∂λX0
This gives
F)F+
1 X λ - λV 2F XV
(
)
}
or
FI ) 0
(8)
Case II: Fixed final time, free final volume.
λX(tf) + V(tf) ) 0 λV(tf) + X(tf) ) 0
}
or
FII ) 0
(9)
Case III: Fixed final volume, free final time.
λX(tf) + V(tf) ) 0 V(tf) - Vf ) 0 H(tf) ) 0
}
or
FIII ) 0
(10)
Case IV: Free final volume, free final time.
λX(tf) + V(tf) ) 0 λV(tf) + X(tf) ) 0 H(tf) ) 0
}
or
FIV ) 0
[
(7)
In this work we analyze the optimum profile for the four cases mentioned earlier. For each of these cases eqs 3(a, b), 4b, and 6(a, b) have to be solved such that they satisfy certain conditions. We now discuss these conditions for each of the four cases. Case I: Fixed final volume fixed final time. The optimum profile for this case is obtained when the following conditions are satisfied
λX(tf) + Vf ) 0 V(tf) - Vf ) 0
[ ]
UIII ) UIV ) [λX0, λV0, tf]t
(11)
The optimal profile is obtained for each case when the respective conditions are satisfied. The optimization problem is a two-point boundary value problem because the conditions on X and V are defined at t ) 0, whereas those on the adjoint variables have to be satisfied at t ) tf. The solution hence has to be obtained iteratively using a technique like the Newton Raphson method. Thus, for cases I and II we iterate on variables λX0 and λV0, whereas for cases III and IV we iterate on λX0, λV0 and tf. The extension of the Newton Raphson method for solving algebraic equations to two-point boundary value problems is called the shooting method. We now discuss this method for case II, the extension to other cases being trivial. Equations 3a and b and 6a and b are integrated from t ) 0 to tf. The initial conditions on X and V are known. For the integration we guess the initial values λX0 and λV0. After the integration, we verify if the condition specified by eq 9 is satisfied. If it is not satisfied, we used the iterative scheme
∂λX ∂λV0 ∂V ∂λV0
(12)
∂λX ∂V + ∂λX0 ∂λX0 JII ) ∂λ ∂X V + ∂λX0 ∂λX0
∂λX ∂V + ∂λV0 ∂λV0 ∂λV ∂X + ∂λV0 ∂λV0
∂λX ∂λX0 ∂V JIII ) ∂λX0 ∂H ∂λX0
∂λX ∂t ∂V ∂t ∂H ∂t
[
[ ]
∂λX ∂V + ∂λX0 ∂λX0 ∂λV ∂X JIV ) + ∂λX0 ∂λX0 ∂H ∂λX0
∂λX ∂λV0 ∂V ∂λV0 ∂H ∂λV0
∂λX ∂V + ∂λV0 ∂λV0 ∂λV ∂X + ∂λV0 ∂λV0 ∂H ∂λV0
]
∂λX ∂V + ∂t ∂t ∂λV ∂X + ∂t ∂t ∂H ∂t
(13)
(14)
]
(15)
The iterative procedure is repeated until we achieve convergence. The routine ‘ode45’ of MATLAB was used for the integrations. In the scheme just presented, the elements of the Jacobian matrix have to be evaluated at the terminal stage (i.e., final time). Hence, the solution technique relies on the accurate evaluation of these derivatives. The equations that govern the evolution of these derivatives are given in Appendixes I and II. Results and Discussion The method of solution just discussed was used to obtain the performance of the fed-batch reactor for all four cases. The kinetic parameters chosen for all the computations reported in this paper are µ ) 5.2 h-1, R ) 1.0 L/g, β ) 25.0 and γ ) 62.9 g/L. The other parameters chosen for the optimization are V0 ) 7 L, S0 ) 0.143 × 10-3 g/L, F ) 4, X0 ) 0.1 g/L, SF ) 500 g/L, Vf ) 10 L, Y ) 0.8, and tf ) 10 h. For this set of parameters, the corresponding optimal flow rate profile for case I is shown in Figure 1. The biomass (XfVf) at the end of the cycle is 1.243. The biomass obtained when the substrate is added at an average flow rate F h ) 0.3 L/h is XfVf ) 1.226. The biomass while operated in the batch mode (i.e., all the substrate is added initially to the reactor) is 0.9416. For case I, there is thus an improvement of 1.375% and 24.25%, respectively, over fed-batch operation with a constant flow rate and batch operation. The effect of the parameter F on the optimal feed rate policy is depicted in Figure 2. For high values of F, the
Ind. Eng. Chem. Res., Vol. 37, No. 11, 1998 4317 Table 1. Comparison of Performance (Xf, Vf) for Different Cases Batch, Semi-Batch (Constant F), and Optimal Fed-Batcha parameter case
t0
tf
V0
Vf
Xf Vf fed-batch
Xf Vf semi-batch with constant F
Xf Vf batch
fixed final time, fixed final volume fixed final time, free final volume free final time, fixed final volume free final time, free final volume
0 0 0 0
10.0 10.0 11.3375 32.2009
7.0 7.0 7.0 7.0
10.0 9.8417 10.0 15.0809
1.2430 1.2636 1.3518 2.5171
1.2259 1.2437 1.3213 2.1317
0.9416 0.9508 0.9796 1.2393
a
X0 ) 0.1, S0 ) 0.143 × 10-3, F ) 4.0, SF ) 500, F h ) 0.3.
Table 2. Effect of G on the Biomass Production (Xf, Vf) in Fed-Batch Mode Comparison with Semi-Batch (Constant F) and Batcha parameter case
r
t0
tf
V0
Vf
Xf Vf fed-batch
Xf Vf semi-batch with constant F
Xf Vf batch
fixed final time, fixed final volume fixed final time, free final volume free final time, fixed final volume free final time, free final volume
0.8 2.0 3.0 3.8
0 0 0 0
10.0 10.0 12.2464 28.2330
7.0 7.0 7.0 7.0
10.0 9.5798 10.0 13.6668
1.4081 1.3393 1.4554 2.5578
1.2259 1.2769 1.3904 2.0465
0.9416 0.9686 1.0064 1.2078
a
X0 ) 0.1, S0 ) 0.143 × 10-3, SF ) 500, F h ) 0.3.
Figure 1. Optimal flow rate profile for case I for F ) 4.0 for biomass production (parameters as in text).
Figure 2. Effect of F on the optimal flow rate for maximizing biomass production for case I.
deviation from the average flow rate has more weightage in the objective function eq 4b. Hence, this has to be minimized and so the optimal flow rate profile is as close to F h as possible. Consequently, the improvement over the constant flow rate is marginal. This yields XfVf ) 1.243 for F ) 4.0. When we decrease the parameter F, the optimal profile exhibits a significant deviation from F h . Now the weightage of the cost of the deviation from F h has been decreased. So, the optimal profile can deviate from F h and there is a substantial improvement in the total amount of biomass produced (for F ) 0.8, XfVf ) 1.4081). Quantitatively, we now have a 15% improvement over the constant flow rate case.
The results for the other four cases of operation are summarized in Table 1. In this table we have used the same weighting factor F as well as the same feed and initial substrate concentrations for all four cases. The improvements in biomass produced for case IV are 50.76% and 15.31%, respectively, over batch and a fedbatch operation with a constant flow rate. This is the total optimization problem because here we constrain neither the final time nor the final volume. The improvements in reactor performance for cases II and III, where we constrain either the final volume or the final time, lie in between those for case I, where both are constrained, and case IV, where there are no constraints. This table quantitatively verifies the results, which we expect intuitively; that is, the unconstrained optimization problem possesses a much higher optima than the most constrained problem. As discussed earlier, the improvement in reactor performance as measured by the amount of biomass produced at the end of the cycle is more pronounced as we decrease F. However as we decrease F, the problem becomes singular and the convergence of the numerical scheme is not assured. In Table 2 we present numerical results of our simulations of the optimization problems for the four cases for a different set of parameter F. In this table, the values of F are chosen is such that they are the lowest values for which the numerical scheme converges with ease. For values of F lower than those indicated, the convergence is difficult. We have shown results for the lowest value because this is the one that gives the most significant improvement over the fedbatch operation with constant F h and the batch operation. The effect of the parameter F on the optimal performance can be seen by comparing the results presented in Table 2 and Table 1. As F is decreased, the total biomass produced at the end of the reactor cycle increases for each case. The corresponding optimal flow rate profile for each of these cases is depicted in Figure 3 (a-d). The values of the parameter F for which the simulations were performed are also indicated in Table 2. In the nonrepeated fed-batch mode of operation, the optimal profile is dependent on the initial conditions. The kinetic expression is characterized by the substrate concentration S* at which the reaction rate is a maxi-
4318 Ind. Eng. Chem. Res., Vol. 37, No. 11, 1998
Figure 3. Optimal flow rate profiles (parameters as in Table 2): (a) case I; (b) case II; (c) case III; (d) case IV.
mum. Two other characteristic quantities quantifying the substrate concentration are S0, the initial concentration, and SF, the feed concentration. Cazzador (1988) determined the effect of initial conditions on the optimal flow rate profile analytically. His technique, however, cannot be easily extended to more complex systems easily. We now discuss this effect for different ordering of S0, S*, and SF. S* < SF < S0. The optimal flow rate profile for this case is shown in Figure 4a. The flow rate decreases with time and the initial substrate concentration is very large. The reaction rate at this S0 is low. Starting with a high ‘F’, we reduce the substrate concentration from S0 to as close to S* as possible, which increases the reaction rate and results in an optimal performance. S* < S0 < SF. Here the optimal flow rate profile increases with time (Figure 4b). The reaction rate at S0 is lower than the maximum. Addition of substrate at SF would serve only to increase S and thereby decrease the reaction rate further. The addition of the substrate, hence, has to be low initially, as the addition of high F(t) initially would increase S0 to a value away from S*, thereby lowering the reaction rate. It is necessary to add substrate as we have to satisfy the constraint on the final volume. S0 < S* < SF. Here again, the flow rate profile decreases with time. The addition of a significant amount of substrate (high SF) initially increases S0 close to S*, which increases the reaction rate. As time progresses, the flow rate is decreased so as to maintain SF close to S* (Figure 4c). S0 < S* < SF. This is the same as the ordering just discussed with SF being very high. Here, however, because SF is significantly high, the flow rate is initially
high. This raises S0 to S* the flow rate is then decreased to prevent SF from increasing beyond S*, which would decrease the reaction rate. The gradual increase in the flow rate in the second half ensures that the substrate concentration is around S* (Figure 4d). SF < S* < S0. Here, the initial concentration is high so the reaction rate is low. A significant amount of SF is added initially, thereby reducing S0 to close to S*, which generates the decreasing profile of F(t) (Figure 4e). S0 < SF < S*. Here, S0 is low and by adding SF at a high flow rate we increase S0 closer to S*, thereby increasing the reaction rate. The flow rate then decreases with time (Figure 4f). SF < S0 < S*. Here, the initial concentration is lower than S*. The optimal feed rate profile increases with time. Here, the initial concentration is closer to S* than SF. Adding a high amount of fresh feed would only lower S0 and the reaction rate, so it is necessary to start with a low flow rate. The flow rate increases with time, thereby filling up the reactor (Figure 4g). The quantitative results of a typical set of calculations are tabulated in Table 3. Simulations were performed for other cases to study the effect of initial conditions. These results also show a similar trend. In most cases, the biomass increases monotonically with time, for some cases, there is an increase followed by a decrease, and for some there is a decrease followed by an increase under optimal conditions. A typical variation of the biomass concentration with time is shown in Figure 5. Conclusions In this work we have obtained optimal flow rate profiles for a singular problem by converting it to a
Ind. Eng. Chem. Res., Vol. 37, No. 11, 1998 4319
Figure 4. Optimal flow rate profiles (parameters as in Table 3): (a) S* < SF < S0; (b) S* < S0 < SF; (c) S0 < S* < SF; (d) S0 < S* < SF; (e) SF < S* < S0; (f) S0 < SF < S*; (g) SF < S0 < S*.
nonsingular problem. This conversion is achieved by modifying the objective function by including a factor representative of the cost for implementing the optimal
profile. As we decrease the weighting parameter F, the importance on the variation of F is reduced. Then, a more significant variation in F is permissible and there
4320 Ind. Eng. Chem. Res., Vol. 37, No. 11, 1998 Table 3. Effect of Initial Conditions on the Flow Rate Profile and Biomass Production for Fixed Final Time, Fixed Final Volume parameter case 1. 2. 3. 4. 5. 6. 7.
S* < SF < So S* < So < SF S0 < S* < SF S0 < S* < SF SF < S* < So S0 < SF < S* SF < So < S* a
F
So
Sf
0.05 15.0 10.0 0.05 10.0 15.0 -3 0.8 0.143 × 10 25.0 4.0 0.143 × 10-3 500 0.8 25.0 0.143 × 10-3 0.3 0.143 × 10-3 5.0 0.3 5.0 0.143 × 10-3
XfVf fed-batch 2.3140 2.4454 2.1429 1.2430 2.0394 1.2113 2.2492
Xo ) 0.1, tf ) 10.0, V0 ) 7.0, Vf ) 10.0, S* ) 7.9310.
is not based on any laboratory experiments. It has implications in the design of optimal systems where cost implications have to be considered for the maintenance of the optimal profile of the manipulated variable the. The techniques proposed by Modak (1986) could then be used to numerically solve the problem elegantly. In the future, we plan to implement this strategy to solve more complex problems that are multidimensional. Appendix I Equation 7 is substituted back in eqs 3a, 3b, 5, 6a, and 6b. The sensitivity is determined with respect to λX0 and λV0 by differentiating the evolution equations with respect to λX0 and λV0, for cases I and II. This yields
[
]
∂µ(g(X,V)) F ∂X ∂X˙ ) µ(g(X,V)) + X + ∂λX0 ∂X V ∂λX0 ∂µ(g(X,V)) FX ∂V X ∂F X (1a) + 2 ∂V V ∂λX0 V ∂λX0
{
[
}
]
∂µ(g(X,V)) F ∂X ∂X˙ ) µ(g(X,V)) + X + ∂λV0 ∂X V ∂λV0 ∂µ(g(X,V)) FX ∂V X ∂F X (1b) + 2 ∂V ∂λ V ∂λV0 V V0
{
Figure 5. Variation of biomass concentration with time for the flow rate profile in Figure 3a.
is a significant improvement in biomass produced. The convergence of the numerical scheme is ensured for large values of F. As we lower the parameter F, however, the optimization program approaches the singular limit. For low F, the optimization program does not converge easily. The cutoff value of F, above which convergence occurs, depends on the case being studied. This cutoff value is low for case I and high for case IV. This result is to be expected because the case IV is the completely optimal problem because there are no constraints on the final volume and final time. The optimal profile generated using our approach exhibits a continuous variation in F(t). This result is to be expected because now we have a nonsingular problem. The optimal profiles generated can be used to obtain qualitative insight into the sequence of control actions of the singular problem. Thus, for the case being present in Figure 3a, the optimal profile for the singular problem would consist of a time interval where F(t) is Fmax followed by a period where it is Fmin. A singular arc would then follow this in turn. The advantage of the technique proposed is based on the fact that the iterative scheme integrates only a system of ordinary differential equations up to the final time. This happens because the problem is nonsingular and the integrations are performed using F as defined in eq 7. This situation is to be contrasted with the integrations of the singular problem where it is necessary to check if some conditions are satisfied and depending on this integrate different systems of equations. The technique proposed in this work can be used to obtain qualitative insight in the different sequences of control actions for a singular problem. The work discussed here is primarily of a theoretical nature and
}
[
∂F ∂V˙ ) ∂λX0 ∂λX0
(1c)
∂V˙ ∂F ) ∂λV0 ∂λV0
(1d)
] ]
∂λ˙ X ∂2µ(g(X,V)) ∂X ∂µ(g(X,V)) +X ) - λX 2 ∂λX0 ∂X ∂λX0 ∂X2
[
∂µ(g(X,V)) ∂2µ(g(X,V)) F ∂V +X + 2 ∂V ∂V ∂X V ∂λX0 λX ∂F ∂µ(g(X,V)) F ∂λX µ(g(X,V)) + X + (1e) ∂X V ∂λV0 V ∂λX0 λX
[
]
[
] ]
∂λ˙ X ∂2µ(g(X,V)) ∂X ∂µ(g(X,V)) +X ) - λX 2 ∂λV0 ∂V ∂λV0 ∂X2
[
2
∂µ(g(X,V)) ∂ µ(g(X,V)) F ∂V +X + 2 ∂V ∂V ∂X V ∂λV0 λX ∂λX ∂µ(g(X,V)) F ∂λX µ(g(X,V)) + X + (1f) ∂X V ∂λV0 V ∂λV0
[
λX
[
]
]
∂λ˙ V ∂µ(g(X,V)) ∂2µ(g(X,V)) F ∂X +X + 2 ) - λX ∂λX0 ∂V ∂X∂V V ∂λV0
[ [
2
∂ µ(g(X,V))
[
]
2FX ∂V - 3 ∂V2 V ∂λX0 λXX ∂F X∂µ(g(X,V)) FX ∂λX + 2 - 2 (1g) ∂V V ∂λX0 V ∂λX0
λX X
]
]
∂λ˙ V ∂µ(g(X,V)) ∂2µ(g(X,V)) F ∂X +X + 2 ) - λX ∂λV0 ∂V ∂X∂V V ∂λV0
[ [
2
∂ µ(g(X,V))
]
2FX ∂V ∂V V3 ∂λV0 λXX ∂F ∂µ(g(X,V)) FX ∂λX + 2 X - 2 (1h) ∂V V ∂λV0 V ∂λV0
λX X
2
-
]
Ind. Eng. Chem. Res., Vol. 37, No. 11, 1998 4321
The equations for ∂F/∂λX0 and ∂F/∂λV0 are given by
[( ) [( )
( ) ( )
] ]
∂λV λXX ∂V X ∂λX 1 λX ∂X ∂F ) + 2 ∂λX0 2F V ∂λX0 V ∂λX0 ∂λX0 V ∂λX0
()
(1j)
∂λV λXX ∂V X ∂λX ∂F 1 λX ∂X ) + (1k) 2 ∂λV0 2F V ∂λV0 V ∂λV0 ∂λV0 V ∂λX0
()
Appendix II The sensitivity equations required for case 3 and 4, in addition to the equations given in Appendix I, are:
[
]
[
FX ∂λX ∂H ) Xµ(g(X,V)) + λX µ(g(X,V)) + ∂λX0 V ∂λX0 ∂µ(g(X,V)) FX ∂V ∂µ(g(X,V)) F ∂X + λX X + X + 2 ∂X V ∂λX0 ∂V V ∂λX0 λXX ∂λV ∂F + (1l) F + λV + 2F(F - F) ∂λX0 V ∂λX0
[
]
[
]
[
]
[
]
∂H FX ∂λX ) Xµ(g(X,V)) + λX µ(g(X,V)) + ∂λV0 V ∂λV0 ∂µ(g(X,V)) FX ∂V ∂µ(g(X,V)) F ∂X + λX X + X + 2 ∂X V ∂λV0 ∂V V ∂λV0 ∂λV λXX ∂F F + λV + 2F(F - F) + (1m) ∂λV0 V ∂λV0
[
]
[
]
[
[
]
]
∂H FX ∂λX ) Xµ(g(X,V)) + λX µ(g(X,V)) + ∂t V ∂t ∂µ(g(X,V)) F ∂X ∂µ(g(X,V)) FX ∂V X + λX X + 2 + ∂X V ∂t ∂V V ∂t ∂λV λXX ∂F F + + λV + 2F(F - F) (1n) ∂t V ∂t
]
[
[
Nomenclature F: flow rate (L/h) Fmax: maximum flow rate (L/h) Fmin: minimum flow rate (L/h) H: Hamiltonian So: initial substrate concentration (g/L) SF: feed substrate concentration (g/L) to: initial time (h) tf: finial time (h) Vo: initial volume (L) Vf: final volume (L) Xo: initial biomass concentration (g/L) Xf: final biomass concentration (g/L) Y: yield of the reaction Greek Letters µ: kinetic expression
]
]
µmax: maximum specific growth rate h-1 R, β, γ: parametric constants in the kinetic expression λx, λv: adjoint variables in X,V F: weightage of cost factor Subscripts F: feed f: final o: initial max: maximum min: minimum
Literature Cited Biegler, L. T. Solution of Dynamic Optimisation, Problems by Successive Quadratic Programming and Orthogonal Collocation. Computers Chem. Eng. 1984, 8, (3/4), 243-248. Bonte, P.; Modak, J. M.; Tayeb, Y. J.; Lim, H. C. Computational Algorithms for Optimal Fed Rates for a Class of Fed Batch Fermentation: Numerical Results for Penicillin and Cell Mass Production. Biotechnol. Bioeng. 1996, 28, 1408-1420. Cazzador, L. On the Optimal Control of Fed-Batch Reactors with Substrate-Inhibited Kinetics. Biotechnol. Bioeng. 1988, 31, 670674. Cuthrell, J. E.; Biegler, L. T. Simultaneous Optimization and Solution Method for Batch Reactor Control Profiles. Computers Chem. Eng. 1989, 13, (1), 49-63. Isermann, R. Digital Control Systems; Springer Verlag: Berlin, 1981. Lee, J.; Parulekar, S. H. Enhanced Production of R-Amylase in Fed Batch Cultures of Bacillus Subtilis TN106. Biotechnol. Bioeng. 1993, 42, 1142-1150. Levin, K. L. Maximising the Product Distribution in Batch Reactors: Reactions in Parallel. Chem. Eng. Sci. 1992, 47, 1751-1760. Leis, J. R.; Kramer, M. A. Sensitivity Analysis of System of Differential and Algebraic Equations. Computers Chem. Eng. 1985, 9, (1), 93-96. Morrison, K. R.; Sargent, R. W. H. Optimisation of Multistage Processes Described by Differential-Algebraic Equations. Lecture Notes in Mathematics 1984, 12, 86-102. Modak, J. M.; Lim, H. C.; Jayeb, Y. K. General Characteristics of Optimal Feed Rate Profiles for Various Fed Batch Fermentation Processes. Biotechnol. Bioeng. 1986, 28, 1396-1407. Shukla, P. K.; Pushpavanam, S. Optimisation of Biochemical Reactors: An Analysis of Different Approximations of Fed Batch Operation. Chem. Eng. Sci. 1998, 53, 341-352. Vassiliadis, V. S.; Sargent, R. W. H.; Pantelides, C. C. Solution of a Class of Multistage Dynamic Optimisation Problem 1: Problems Without Path Constraints. Ind. Eng. Chem. Res. 1984, 33, 2111-2122. Weigand, W. A. Maximum Cell Productivity by Repeated FedBatch Culture Yield Case. Biotechnol. Bioeng. 1981, 23, 249266. Westage, J. P.; Curtis, R. W.; Emergy, H. A.; Hasegawa, M. P.; Heinstein, P. F. Approximation of Continuous Growth of Cephalotaxus Harringtonia Plant Cell Cultures Using Fed Batch Concentrations. Biotechnol. Bioeng. 1991, 38, 241-246.
Received for review August 19, 1997 Revised manuscript received July 14, 1998 Accepted July 15, 1998 IE970557W
