Ind. Eng. Chem. Res. 2004, 43, 4211-4219
4211
PROCESS DESIGN AND CONTROL Batch Fermentation Networks Model for Optimal Synthesis, Design, and Operation Gabriela Corsano,† Pı´o A. Aguirre,*,† Oscar A. Iribarren,‡ and Jorge M. Montagna‡ INGAR, Instituto de Desarrollo y Disen˜ o, Avellaneda 3657, (3000) Santa Fe, Argentina
This paper addresses an integral optimization of fermentation processes. The behavior of the fermentors is described by a set of algebraic and differential equations written as finite-difference equations in an equation-oriented environment. Unconventional constraints related to the number of batch items and connections among them, detailed kinetic models and operating costs corresponding to inoculum, and different available substrates are included in the model. The optimal number of units to be used in the process, their optimal operation policy (i.e., connected in series or in parallel working out of phase), as well as the optimal volume and operation of each unit, are determined simultaneously. The model is formulated as a sequence of nonlinear programming (NLP) problems. 1. Introduction In recent years, there has been significant growth in chemical industry in the use of batch fermentors because of the demand for a large number of specialty chemicals. Genetically engineered microorganisms produce complex molecules of higher quality more efficiently than does chemical synthesis. Fermentation isalso useful in the case of less sophisticated products, as it allows their production starting from chemically complex but rather inexpensive raw materials that are byproducts or even wastes from agricultural processes. In general, batch processing is used in manufacturing low-volume, high-value products. Thus, even a moderate increase in product yield can lead to a considerable improvement in profitability. For this reason, it is important to address the modeling and optimization of these batch processes.1 There is abundant literature on batch process synthesis and design with batch stages described by fixed time and size factors as in refs 2-5. The papers by Reklaitis and co-workers2,3 resorted to algorithmic solution procedures, whereas the latter works used the mixed-integer nonlinear programming (MINLP) approach: for batch processes in general in Ravemark and Rippin4 and for processes that specifically include fermentation stages in Montagna et al.5 This approach allows for the optimization of plant decision variables, including batch sizes and operating times of semicontinuous items, as well as the structure of the plant, i.e., the number of units in parallel and the provision of storage tanks between stages. However, the use of constant time and size factors requires that the units’ process decision variables, e.g., * To whom correspondence should be addressed. Tel.: 54-342-4534451. Fax: 54-342-4553439. E-mail: paguir@ ceride.gov.ar. † Universidad Nacional del Litoral. ‡ Universidad Tecnolo ´ gica Nacional.
reaction extents, be fixed, thus preventing better solutions for these problems. A first level of detailed description of the units’ performance depending on these process variables consists of using algebraic models. Such an approach was first proposed by Salomone and Iribarren6 and applied to fermentation processes by Pinto et al.7 This approach allows for the simultaneous optimization of the units’ process variables and the plant decision variables. The process performance models are additional algebraic equations describing the time and size factors as functions of the units’ process variables, so that the global problem formulation is still an MINLP. A more detailed description of the performance of batch stages requires that they be modeled with differential equations. First Barrera and Evans8 and later Salomone et al.9 proposed that this simultaneous optimization should be approached by integrating the batch plant model with dynamic simulation modules for the batch units. This approach does succeed in performing the optimization, but it entails a great computational effort. The global problem is no longer just an MINLP. It has dynamic simulation blocks interacting with it, and therefore, this requires algorithmic solution procedures. Incidentally, the approach in ref 9 proposes an algorithm whose resolution sequence overcomes the unfeasibility problems reported in ref 8. One way of incorporating the units’ dynamic models without losing the MINLP nature of the global problem is to discretize the differential equations to convert them into algebraic constraints of the program. This was the approach used by Bathia and Biegler10 for simple process examples, and it is used here for a complex fermentation network. One of the main motivations for using a detailed model for the fermentation network was to be able to optimize the consumption of complex substrate mixtures that can be fed at different locations of the network. This issue arises in our case study, which is a fermentation
10.1021/ie030549h CCC: $27.50 © 2004 American Chemical Society Published on Web 06/24/2004
4212 Ind. Eng. Chem. Res., Vol. 43, No. 15, 2004
Figure 1. Fermentation process scheme.
process in which subproducts or even residues of a food process are used as substrates for the fermentations. Another point in which we were interested was the inclusion of the amount of inoculum as an optimization variable. This introduces a structural optimization, which is the number of fermentation stages connected in series: use of more units requires less inoculum to be fed to the initial, smaller prefermentor. Overall, an original single fermentation step could be split into several steps connected either in series or in parallel working out of phase or a combination of these alternatives. In this paper, we consider the optimization of fermentation using highly structured models. The technique proposes a special superstructure that is optimized to solve the fermentor network synthesis problem. A network of fermentors involving a large number of units is postulated, and the optimal subnetwork that optimizes the specified objective function is found. The number of fermentors to be used in this process, their configuration, and their size are obtained. The design variables of the batch stages, namely, size factors, processing times, and operating costs, are computed as functions of the process variables. In this way, the resolution of the synthesis problem is performed simultaneously, taking into account fermentors model equations and their interconnections as NLP problem constraints. Even if structured modeling can be both data- and computer-intensive,11 we show that the problem at hand can be solved efficiently. To our knowledge, there is no previously published work dealing with the optimal design, based on such a complete model, of the batch fermentation network. The remainder of the paper is organized as follows: The next section describes a general fermentor model. In section 3, we examine the superstructure for efficient synthesis problems considering different scheduling scenarios. In section 4, numerical examples for ethanol and torula yeast productions are considered: we present the results and discuss the optimal synthesis, design, and operation of these production processes. Finally, conclusions of this work are outlined in section 5. 2. Fermentor Model The proposed models consists of a set of equations that describe mass balances in the fermentors, fermentors connection constraints, scheduling constraints, and the objective function. Figure 1 shows, as an example, a simple fermentation network that consists of units in series. These fermentors could be fed with molasses diluted with water, mill juices, filters juice, and/or distillery vinasses. Molasses is a byproduct obtained from the sugar plant crystallizer (part of the mother liquor, which is purged to get rid of components different from sugar), whereas the mill and filters juices are intermediate products of the same process. Distillery vinasses or distillery broth is the nondistilled residue of the ethanol process. We limited this study to batch fermentors, so no substrates are added during the course of fermentations.
In the case of biomass-producing fermentors connected in series, the product of one fermentor is the inoculum of the next. In addition, we assume that the amount of inoculum added to the first fermentor of the network (a broth containing biomass prepared in the laboratory) is a decision variable and its cost is included in the overall objective function. The symbols used to describe the process models are defined in the Nomenclature section. 2.1. Mass Balances. The fermentation balances for a fermentor f can be expressed by the following differential equations using Monod’s kinetic model12
Biomass
dXf ) µfXf - υfXf dt
(1)
Substrate
dSf µfXf )dt Yf,x/s
(2)
Dead Biomass
dXdead f ) υfXf dt
(3)
Product
dPf µfXf ) dt Yp/x
(4)
Sf µf ) µmax,f ks + Sf
(5)
where
The yield coefficient, Yf,x/s, which is an efficiency measure for a particular conversion, is in this case the substrate-biomass conversion and has been experimentally found to depend on the quality of the carbohydrate source. In this work, we use the following expression for this coefficient, obtained through a fit to industrial data13
Yf,x/s ) 0.3689x1,f + 0.4092x2,f + 0.2522x3,f + 0.1792x1,fx2,f + 0.3736x1,fx3,f + 0.1123x2,fx3,f (6) where
x1,f + x2,f + x3,f ) 1
(7)
This means that the yield coefficient has a nonlinear dependence on the contributions of the molasses, juices, and vinasses. 2.2. Finite-Difference Equations. We use the trapezoidal method to discretize differential eqs 1-4. Let h g 0 and define the nodes by tn ) t0 + nh, n g 0. The general form of the trapezoidal method to solve dx/dt ) f(t,x) is
h xn+1 ) xn + [f(tn,xn) + f(tn+1,xn+1)] 2 This is an implicit one-step method, and it exhibits a special stability property.14 We solve the differential equations expressed by this method in GAMS15 as a set of constraints of the optimization problem. The performance of this method was tested with the gPROMS simulator16,17 for different productions and different numbers of node points. Figure 2 shows an example of biomass growth for yeast production, obtained with 5 and 30 nodes in the trapezoidal discretization solved in GAMS with the gPROMS simulator.
Ind. Eng. Chem. Res., Vol. 43, No. 15, 2004 4213
Figure 2. Trapezoidal method in GAMS versus gPROMS simulator.
As can be observed, the trapezoidal method exactly matches the simulation solution at the last grid point (solution point), and this occurs in all simulated cases (biomass concentration, substrate concentration, product concentration) with various numbers of grid points. This is an advantage of solving the differential equations inside an optimization problem because, near the optimal solution, the fermentation process variables taper off, and thus gradient elements at these points are small. Therefore, the error at the last point of the trapezoidal method solved in GAMS is negligible, even for as few as five grid points. Equations 8-11 represent the discretization of differential eqs 1-4
Xf,i+1 ) Xf,i +
{[ ( ) ] [ ( ) ] } [ ( ) ( ) ]
hf Sf,i µmax,f - υf Xf,i + 2 Sf,i + ks Sf,i+1 µmax,f - υf Xf,i+1 Sf,i+1 + ks
Pf,i+1 ) Pf,i +
[ (
) (
productivity )
(10)
hf Sf,i Xf,i µ + 2 max,f Sf,i + ks Yx/p Sf,i+1 Xf,i+1 (11) µmax,f Sf,i+1 + ks Yx/p
) ]
∀f ) 1, ..., nfer; ∀i ) 1, ..., nnodes - 1 where nfer represents the number of fermentors and nnodes to the number of nodes of the trapezoidal discretization. Thus Xf,0 is the initial biomass concentration, Xf,nnodes is the final biomass concentration, Sf,0 is the initial substrate concentration, and so on. 2.3. Scheduling Equations. We computed two sets of scheduling constraints: (a) relaxed problem and (b) zero-wait (ZW) transfer. The first set assumes that the batch can be stored so that each piece of equipment can
Pnfer,nnodesVnfer Tnfer
(12)
The ZW policy considers that a processed batch is immediately transferred to the following stage. For this policy, we consider adding parallel units. In this case, the plant cycle time is equal to the maximum stage cycle time, which is computed as the operating time divided by the number of parallel units at this stage, that is
CT ) max (Tf/NPf) f)1,...,nfer
(8)
hf Sf,i Xf,i Sf,i+1 ) Sf,i + -µmax,f 2 Sf,i + ks Yf,x/s Sf,i+1 Xf,i+1 (9) µmax,f Sf,i+1 + ks Yf,x/s hf dead dead ) Xf,i + (υfXf,i + υfXf,i+1) Xf,i+1 2
operate at its own optimal batch size and cycle time, as if we were using storage enough to decouple the stages, for example. Actually, we do not consider intermediate storage to be operationally feasible in our case study examples. Because of biomass degradation during holding times at tanks, the transfer policy of choice for these cases is ZW, whereas for the anaerobic production of ethanol, storage tanks would behave as fermentors as well. The decoupling intermediate storage policy, without computation of a cost for the storage involved, is used as a constraint-free scenario, as presented in Asenjo et al.18 As this scenario removes all batch plant scheduling constraints, we can solve the optimization of process variables unbiased with respect to the plant structure. In particular, in the absence of scheduling constraints, there is no driving force for duplicating units in parallel, so we use this scenario (the superstructure of which is represented by Figure 1) to optimize the number of stages in series, i.e., the tradeoff between the lower total fermentor volume associated with the larger the number of fermentors in series, each one incorporating additional substrate, and the economic penalty of an increased number of units due to economies of scale. For this problem, the production rate constraint is
(13)
Because of the superstructure representation proposed in this work, the definition of discrete variables is not needed, so eq 13 was reformulated. For this policy, the production rate constraint is
productivity )
Pnfer,nnodesVnfer CT
(14)
2.4. Feeding Equations and Fermentor Connection Constraints. Before beginning each fermentation, the units are fed with either molasses, filter juices, vinasses and water, or a blending of these single feeds. Here, we present the feeding equations and fermentor connection constraints for the two mathematical problems that are solved repeatedly in the case studies: the relaxed problem (RP) and the zero-wait (ZW) transfer policy. Volume Equations.
Mf + FJf + DVf + Wf ) FEEDf
∀f ) 1, ..., nfer (15)
For the RP
Vf-1 FEEDf Vf + ) Tf-1 Tf Tf
∀f ) 2, ..., nfer
(16)
4214 Ind. Eng. Chem. Res., Vol. 43, No. 15, 2004
For the ZW policy
Substrate
Vf-1 + FEEDf ) Vf
∀f ) 2, ..., nfer
(17)
For f ) 1, we consider the inoculum volume
FEED1 + Vinoc ) V1
Sf-1,nnodesVf-1
Sf,FEEDFEEDf Sf,0Vf ) Tf Tf S V + S FEED ) S For the ZW Policy f-1,nnodes f-1 f,FEED f f,0Vf
For the RP
Tf-1
+
∀f ) 2, ..., nfer (26)
(18) Product
where
Vinoc ) Winoc +
Inoc Finoc
(19)
Pf-1,nnodesVf-1
Pf,0Vf Tf V ) P P For the ZW Policy f-1,nnodes f-1 f,0Vf
For the RP
Tf-1
)
Mass Balances. Feeding Sugar Contribution Equations.
MfSmolass + FJfSFJ + DVfSDV ) FEEDfSf,FEED For the RP
x1,fVfSf,0 MfSmolass x1,f-1Vf-1Sf-1,nnodes ) + Tf Tf Tf-1 V x x2,fVfSf,0 FJfSFJ 2,f-1 f-1Sf-1,nnodes ) + Tf Tf Tf-1 V x x3,fVfSf,0 DVfSDV 3,f-1 f-1Sf-1,nnodes ) + Tf Tf Tf-1
(20)
}
x1,fVfSf,0 ) MfSmolass + x1,f-1Vf-1Sf-1,nnodes x2,fVfSf,0 ) FJfSFJ + x2,f-1Vf-1Sf-1,nnodes x3,fVfSf,0 ) DVfSDV + x3,f-1Vf-1Sf-1,nnodes
The initial biomass for the first fermentor is determined by
(28)
nfer
}
For f ) 1
min TAC )
RfVfβ + OP ∑ f)1 f
(29)
where nfer is the number of pieces of equipment used in the fermentation process; Rf and βf are annualized cost coefficients according to ref 19; Vf is the fermentor size; and OP represents operating costs, which include the raw material costs (molasses, filter juices, water) and the inoculum cost, taken as the cost for producing them at the plant laboratory. Distillery vinasses have no cost because they are ethanol plant wastes that would otherwise require additional costs for their disposal. 3. Fermentation Network Superstructure
x1,fV1S1,0 ) M1Smolass x2,fV1S1,0 ) FJ1SFJ
(23)
x3,fV1S1,0 ) DV1SDV Fermentor Connection Constraints. Biomass
Xf-1,nnodesVf-1
Xf,0Vf Tf V ) X X For the ZW Policy f-1,nnodes f-1 f,0Vf Tf-1
)
}
∀f ) 2, ..., nfer (24) Dead Biomass dead V Xf-1,n nodes f-1
dead Xf,0 Vf Tf-1 Tf dead dead X V ) X For the ZW Policy f-1,nnodes f-1 f,0 Vf
For the RP
∀f ) 2, ..., nfer (27)
Note that, from eqs 21 and 7, we obtain eq 26 for the RP case and, from eqs 22 and 7, we obtain eq 26 for the ZW case. 2.5. Objective Function. The objective function to be minimized is the total annualized cost (TAC) of the process
∀f ) 2, ..., nfer (22)
For the RP
}
XinocVinoc ) X1,0V1
∀f ) 2, ..., nfer (21) For the ZW policy
}
)
}
∀f ) 2, ..., nfer (25)
We start our optimization approach by performing the optimal synthesis for the unconstrained RP by modeling a superstructure with as many units in series as indicated by an upper bound: this superstructure is represented by Figure 1. The optimal solution eventually eliminates units of the network, driving their size to zero. In this way, we solve the model as an NLP problem instead of an MINLP problem. Some numerical restrictions were applied to avoid numerically undefined terms in both the equations and their first derivatives. To perform the optimal synthesis for the actually adopted ZW transfer policy, we consider a superstructure with all possible options with units in series or in parallel, with the total number of units of each structure ranging from 1 to an upper bound, which we took here to be the optimal number of units obtained in the RP. Again, we solve the model as a collection of NLP problems instead of implementing an MINLP problem. Each considered structure has its own production rate, and the sum of these production rates is the target productivity. Also, the global objective function is the sum of the contributions of each structure. The optimal solution chooses the best structural option by driving to zero both the sizes of all units and the amounts of
Ind. Eng. Chem. Res., Vol. 43, No. 15, 2004 4215 Table 1. Parameter Values for Both Case Studies torula biomass fermentor µmax ks Yx/s Yx/p υ R19 β19 inoculum cost molasses cost filter juices cost water cost production
0.5 h-1 20 g/L eq 6 0.02 h-1 115 550 0.43 10 $/kg 31.9 $/ton 25.78 $/ton 0.05 $/ton 100 kg/h
ethanol biomass fermentor
ethanol fermentor
0.5 h-1 20 g/L eq 6
0.1 h-1 20 g/L 0.124 0.23 0.02 h-1 0.02 h-1 115 550 0.43 10 $/kg 31.9 $/ton 25.78 $/ton 0.05 $/ton 100 kg/h
inoculum and raw materials that are not involved in the optimal structure. This NLP strategy allows for the explicit exploration of all alternatives, because their number is relatively modest. On the other hand, the available MINLP algorithm (DICOPT) can miss the optimal solution through the application of linear cuts to nonconvex regions. Another advantage is that, in the proposed collection of NLP problems, the user can provide physically meaningful initializations, thus increasing the robustness and usefulness of the optimization models. The robust initialization scheme used assumed that all units are active. Available computer codes for solving MINLP do not allow the user to initialize the intermediate NLP problems. This NLP model probably requires some extra work on problem setup with respect to that for an MINLP program involving a unique superstructure. However, the resolution CPU time for solving the NLP is very short, and the resolution is robust. All examples were implemented and solved in GAMS on a Pentium IV, 1.60-GHz PC. The code CONOPT was employed for solving the NLP problems. In the RP model, we consider eqs 8-11 for all fermentation stages involved, plus size eq 16; feeding sugar contribution eqs 20, 21, and 23; fermentor connection constraints 24-27; and production rate constraint 12. In the ZW model, we consider the same constraints as indicated in the previous paragraph for each of the considered structures or alternatives, plus eq 13, which corresponds to this policy. 4. Case Study Problems The fermentation network model outlined above was solved for two products production examples: Torula utilis yeast and ethanol production with Saccharomyces cerevisiae. The fermentation model parameters values are listed in Table 1. In both cases, we consider a fixed production rate of 100 kg/h and a time horizon of 7500 h. The contribution of total reducing sugars (TRSs) is equal to 779 g/L for molasses, 100 g/L for filter juices, and 10 g/L for distillery vinasses. We also required that the initial substrate concentration at any fermentor should have an upper bound of 100 g/L to avoid substratum inhibition. 4.1. Torula Yeast Fermentation Process. Yeast production using sugar byproducts is an attractive alternative, mainly if we consider that torula yeast is a source of protein for animal species and is frequently used as protein supplement in their diet. Fermentation is the main operation of this process. In this case, the
Table 2. Optimal Solution: Torula Yeast RP Case Tf (h) fer 1 0.003 fer 2 4.572 fer 3 0.003 fer 4 9.508 fer 5 11.423
Vf (m3)
Xinitial Xfinal (g/L) (g/L)
sugar substrates
Sfinal (g/L)
Yx/s
10-5 39.52 39.54 9.22 0.369 0.008 10.67 47.92 M + FJ + DV 7.22 0.428 10-5 42.16 42.24 99.75 0.374 0.563 1.35 35.85 M + DV 4.91 0.388 29.08 0.83 36.02 M + DV 1.11 0.388
product is the biomass itself (active and dead), so eq 4 is not considered. For this production, we verified that five was a good upper bound on the number of units to be used for the RP. Therefore, the superstructure was formulated with five units, i.e., it could be schematized as in Figure 1, with n ) 5. The optimal solution for the RP is reported in Table 2. It corresponds to the use of three units in series. It is worth noting that the volumes of the first and third units are at their lower bounds, with the ratio of volume to processing time so as to satisfy the production rate constraint. These nonzero solutions at the lower bound imposed on the NLP formulation correspond to the zero solutions of an MINLP formulation, i.e., these units do not exist physically. The optimal solution of the unconstrained RP employed increasing dilution ratios and operating times for the consecutive stages, with a minimum final concentration of substrate at the last stage, i.e., conversion was maximized only at the last stage. The use of alternative fermentation feeds (molasses, filter juices, and vinasses) is also a process optimization variable. These feeds are blended so that the total initial substrate concentration in each fermentor is equal to 100 g/L (upper bound of this process variable). Variable feed blending allows for better use of the units, considering different available feeds with their costs, sugar concentrations, etc. The RP optimal solution employed the concentrated molasses as the main sugar source, diluted with inexpensive distillery vinasses at all stages. The more expensive filter juice is also used but only at the first stage where a small quantity appreciably increased the biomass-substrate yield. For the ZW transfer policy, taking into account that the RP optimal solution includes three units in series, we posed the following alternative structures (Figure 3): option 1, one unit; option 2, two units in series; option 3, three units in series; option 4, two units in parallel at a first stage and one at a second stage; and option 5, one unit at a first stage and two units in parallel at a second stage. Each of these options was modeled, and all of them were embedded into a superstructure model of NLP type. Observe that our underlying assumption is that the optimal design will have a number of units that is less than or equal to the number of units at the optimal solution of the relaxed problem. The optimal solution of this model uses two units in series (option 2 in Figure 3), and the values of the most relevant variables of this option are reported in Table 3. It is worth noting that the optimal solution employs equal operating times for both stages, which avoids idle time in either of them, and a smaller number of units with the purpose of minimizing equipment cost. No filter juice feed was used in this case at the first stage (which is larger than the first stage of the RP optimal solution), at the expense of a smaller yield.
4216 Ind. Eng. Chem. Res., Vol. 43, No. 15, 2004 Table 4. Characteristics of Torula Yeast Superstructure Models model
total annual cost ($)
number of variables
number of constraints
CPU time (s)
255,204
403
377
15.6
256,334
404
382
16.8
RP optimal synthesis ZW optimal synthesis
Table 5. Optimal Solution for Torula Yeast ZW Superstructure Model Varying β cycle time configuration (h)
β 0.6 0.43 0.35 a
3 in series 2 in series 1 fermentor
unit size(s) (m3)
total volume (m3)
unit costa ($)
6.80 0.217, 1.985, 15.87 18.072 850,000 11.75 0.623, 29.99 30.613 593,000 20.59 53.16 53.160 464,000
Unit cost ) R∑jVjβj.
Table 6. Optimal Solutions for Torula Yeast ZW Superstructure Model with Varying µmax
µmax configuration 1 0.5 0.25 a
Figure 3. Torula yeast fermentation process scheme for the ZW model.
Figure 4. Biomass and substrate concentrations for torula yeast ZW case. Table 3. Optimal Solution: Torula Yeast ZW Case Tf (h)
Vf (m3)
fer 1 11.754 0.623 fer 2 11.754 29.99
Xinitial (g/L)
Xfinal (g/L)
0.573 0.733
35.31 36.21
sugar Sfinal substrates (g/L) M + DV M + DV
4.14 1.09
Yx/s 0.388 0.388
Figure 4 shows the optimal biomass and substrate concentrations for the optimal ZW transfer policy solution. Although the dilution ratios at the different stages are now more similar than in the RP solution, the final
2 in series 2 in series 2 in series
cycle time (h)
unit sizes (m3)
total inoculum volume volume (m3) (m3)
5.99 0.264, 15.03 15.294 11.754 0.623, 29.99 30.613 23.51 1.557, 59.94 61.497
0.003 0.009 0.026
unit costa ($) 435,000 593,000 811,000
Unit cost ) R∑jVjβj.
concentration of substrate is again minimal at the last stage. Table 4 reports the costs, numbers of variables and constraints, and CPU times required for the above examples. We used this same process model to carry out a sensitivity analysis, varying the values of parameters away from those held for torula. This allows the parameters with the most significant influence on the process to be identified. These parameters are “critical” in the sense that they produce changes in the process synthesis. Following, we briefly describe the effect of each one of these variations, retaining an upper bound of three units. Equipment Cost Exponent. It is worth noting that, for an exponent equal to 1 (linear cost function), the optimal solution tends to use all of the units of the superstructure, because the optimization problem turns itself into minimization of the total fermentation volume. As the exponent is reduced, the optimal solution tends to use a smaller number of units because a small cost exponent penalizes duplication of units. Table 5 lists the optimal solutions obtained for different values of the cost exponent β. As the cost exponent decreases, the number of units decreases even at the expense of a larger total fermentation volume. Reaction Kinetics. The value of µmax in the torula model was equal to 0.5 h-1. If we consider smaller values, the lower substrate consumption rate implies that longer fermentation times are needed to obtain the productivity target. Table 6 reports different optimal solutions for different kinetics to produce 100 kg/h of yeast. As a result of the very sluggish kinetics, the problem is dominated by equipment cost. Table 6 shows that, in this case, the optimal solution resorts to increasing operating costs, e.g., the amount of inoculum consumed, to reduce the size of the fermentors. Inoculum Cost. If this cost is increased, the optimal solution resorts to placing more units in series so that the first stage (which is the only one fed with inoculum
Ind. Eng. Chem. Res., Vol. 43, No. 15, 2004 4217
Figure 5. Ethanol fermentation process scheme for the RP model. Table 7. Optimal Solutions for Torula Yeast ZW Superstructure Model with Varying Inoculum Cost inoculum cost ($/kg) configuration 1 10 100
1 fermentor 2 in series 3 in series
cycle time (h)
unit size(s) (m3)
11.782 30.105 11.754 0.623, 29.99 11.517 0.011, 0.66, 29.34
inoculum volume (m3)
TAC ($)
0.546 0.009 1.7-4
242,550 256,334 257,175
Table 8. Optimal Solution: Ethanol RP Case Tf (h) biomass fer 1 14.19 biomass fer 2 0.003 ethanol fer 1 24.23 ethanol fer 2 0.003
Vf (m3)
Xinitial (g/L)
Xfinal (g/L)
3.196 0.219 35.023 10-5 35.023 35.03 45.509 4.109 12.012 10-5 12.012 12.011
sugar Sfinal substrates (g/L) M + DV M + DV -
3.988 3.988 3.342 3.338
prepared at the laboratory) is smaller, as shown in Table 7. Then, the final biomass of each unit becomes the inoculum of the next one. 4.2. Ethanol Fermentation Process. The ethanol fermentation process consists of two distinct types of fermentors: the first are of aerated type for biomass production, whereas the latter are anaerobic for ethanol production (eqs 1-4). Therefore, the first stage is similar to that of the torula yeast fermentation process, whereas the second stage has different parameters for the maximum growth rate and yield coefficient (Table 1). After an exploratory study on this fermentation process, we found that the appropriate upper bounds for the RP were two biomass fermentors and two ethanol fermentors, as in the superstructure with up to two units per stage, which is schematized in Figure 5. The RP optimal solution consists of just one biomass fermentation unit and one ethanol fermentation unit. Some of the optimal values for the process variables are listed in Table 8. As can be noted, the second biomass and ethanol fermentors operate for negligibly short times and their sizes are at their lower bounds: this means that they do not exist physically. The alternative structures considered for the ZW transfer policy were (Figure 6): option 1, one biomass fermentor and one ethanol fermentor; option 2, two biomass fermentors in series and one ethanol fermentor; option 3, one biomass fermentor and two ethanol fermentors in series; option 4, two biomass fermentors in parallel and one ethanol fermentor; and option 5, one biomass fermentor and two ethanol fermentors in parallel. Each option was modeled implementing a superstructure model of the NLP type. From the superstructure optimization, we found that the best configuration corresponds to the use of just two units: one biomass fermentor and one ethanol fermentor in series, and again, this happens because of the large penalty for duplication implicit in the cost exponent of 0.43. A sensitivity analysis of the ethanol process was also carried out. First, we present the results of the superstructure optimization with equipment costs being reduced by an order of magnitude i.e., with R ) 11 550, and second, with the cost exponent being increased to β ) 0.6. In both cases, the other parameters are kept constant at their values of Table 1.
Figure 6. Ethanol fermentation process scheme for the ZW model.
Figure 7. Gantt chart for the optimal solution in the ethanol ZW case.
For R ) 11 550, equipment costs are low, so operating costs dominate, and the optimal configuration corresponds to one biomass fermentation unit and two ethanol fermentation units in parallel (option 5 in Figure 6). Considering that the ethanol fermentation time is almost twice the biomass fermentation time in the RP case, this option sounds reasonable: Figure 7 shows the Gantt chart for this optimal solution where the ethanol fermentation time is exactly twice the biomass fermentation time, which allows process variables to attain values close to those they had in the RP optimal solution. A different structure is obtained when we change the size exponent of the cost function. If we take β ) 0.6 with the rest of the parameters corresponding to the original data, the optimal configuration consists of the use of one biomass fermentor and two ethanol fermentors in series (option 3, Figure 6). This happens because
4218 Ind. Eng. Chem. Res., Vol. 43, No. 15, 2004 Table 9. Optimal Solution: Ethanol ZW Case for β ) 0.6
biomass fer ethanol fer 1 ethanol fer 2
Tf (h)
Vf (m3)
Xinitial (g/L)
Xfinal (g/L)
sugar substrates
Pinitial (g/L)
Pfinal (g/L)
Sfinal (g/L)
10.62 10.62 10.62
3.06 8.51 12.40
0.765 12.29 13.83
34.14 20.16 21.76
M + DV M + DV M + DV
0 34.02
49.59 85.61
8.36 8.004 4.30
Table 10. Characteristics of Ethanol Superstructure Models model RP optimal synthesis ZW optimal synthesis ZW optimal synthesisa ZW optimal synthesisb a
total annual number of number of CPU time cost ($) variables constraints (s) 285,900
380
358
6.1
289,875
1157
1090
50.7
73,890
1157
1090
34.7
395,175
1157
1090
24.9
For R ) 11 550, β ) 0.43. b For R ) 115 500, β ) 0.6.
when the exponent increases, it is more economical to use more units with a smaller size than to use fewer units with a larger size. Table 9 lists some of the optimal values of the variables, and Table 10 presents the computational characteristics of the ethanol superstructure model. 5. Conclusions This paper presents a procedure to simultaneously obtain the optimal synthesis, design, and operation of fermentation network processes. In this approach, a general superstructure is proposed and then rigorously modeled as a collection of NLP problems. The formulation involves integration between fermentors as well as variable feeds, sizes, and operating times. The optimal number of units and the best operating policy are determined simultaneously with the optimal values of the process variables. This approach is novel in terms of the decision variables involved, the integration of algebraic and differential equations, and the solution of an NLP problem even when discrete decisions are involved. In particular, we are not aware of previous works dealing with the tradeoff between duplication in series vs duplication in parallel. Numerical examples for different fermentation processes were presented. Small solution times were required, even when highly nonlinear and nonconvex formulations were solved. A sensitivity analysis was carried out to find the critical model parameters, and results of these process parameter variations were presented. In the case of low equipment costs, the optimal solution resorts to units in parallel working out of phase. This is the processing policy that allows process variables to attain values that are closer to those they had at the RP optimal solution, which has no scheduling constraints. In other words, if equipment cost is not relevant, the optimal solution priority is process efficiency in terms of raw materials and inoculum. Even if this is not the case for the considered base cases of the ethanol and torula processes, it will certainly be so in the case of fermentations involving very high valued products. Otherwise, in the case of high equipment costs but with a size exponent closer to unity, the optimal solution resorts to duplication in series, which is the policy that
minimizes the total volume of the units. This would be the case for processes that use sophisticated fermentors of a modular type, as in the case of membrane bioreactors. Raw materials (molasses, filter juices, vinasses, and water) have different costs and substrate concentrations, with the optimal solution choosing the best blend. There exists a tradeoff among carbohydrate sources, concerning the costs of their sugar contributions and their equipment size requirements. In our particular case, optimal solutions used only molasses diluted with vinasses; filter juice and process water were not selected as feeds. Our results show that, depending on the values of the model parameters, different proposed process alternatives can be optimal. This emphasizes the convenience of considering such alternatives in the design. Acknowledgment The authors acknowledge financial support provided by CONICET (Consejo Nacional de Investigaciones Cientı´ficas y Te´cnicas) under Grant PIP 2706, Agencia Nacional de Promocio´n Cientı´fica y Tecnolo´gica under Grant PICT 8752, and UNL (Universidad Nacional del Litoral). Nomenclature DVf ) volume of distillery vinasses added to fermentor f (m3) FEEDf ) feed for fermentor f (m3) FJf ) volume of filter juices added to fermentor f (m3) hf ) step size in trapezoidal method used for fermentor f Inoc ) mass of inoculum (kg) ks ) substrate saturation constant (kg m-3) Mf ) volume of molasses added to fermentor f (m3) nfer ) number of fermentors nnodes ) number of grid points NPf ) number of parallel units in fermentor f Pf ) product concentration in fermentor f (kg m-3) Pf,i ) product concentration in fermentor f at grid point i (kg m-3) SFJ ) total concentration of reducing sugars in filter juice (kg m-3) Smolass ) total concentration of reducing sugars in molasses (kg m-3) Sf ) substrate concentration in fermentor f (kg m-3) Sf,i ) substrate concentration in fermentor f at grid point i (kg m-3) Sfinal ) final substrate concentration (kg m-3) SDV ) total reducing sugars in distillery vinasses (kg m-3) tn ) final time in the trapezoidal method t0 ) initial time in the trapezoidal method Tf ) final time in fermentor f (h) Vinoc ) inoculum size (m3) Vf ) unit size of fermentor f (m3) Wf ) volume of water added to fermentor f (m3) Winoc ) volume of water used for inoculums (m3) x1,f ) percentage contribution of (TRSs) of molasses in fermentor f x2,f ) percentage contribution of (TRSs) of filter juice in fermentor f
Ind. Eng. Chem. Res., Vol. 43, No. 15, 2004 4219 x3,f ) percentage contribution of (TRSs) of distillery vinasses in fermentor f Xf ) biomass concentration in fermentor f (kg m-3) Xf,i ) biomass concentration in fermentor f at grid point i (kg m-3) Xfinal ) final biomass concentration (kg m-3) Xinitial ) initial biomass concentration (kg m-3) Xdead ) inactive biomass concentration in fermentor f (kg f m-3) dead ) inactive biomass concentration in fermentor f at Xf,i grid point i (kg m-3) Yx/p ) product yield coefficient Yf,x/s ) biomass yield coefficient in fermentor f Rf ) cost coefficient for fermentor f βf ) cost exponent for fermentor f µf ) specific growth rate of biomass in fermentor f (h-1 ) µmax,f ) maximum specific growth rate of biomass in fermentor f (h-1 ) Finoc ) inoculum density (kg m-3) υf ) biomass death rate in fermentor f (h-1 )
Literature Cited (1) Rahman, S.; Palanki, S. State Feedback Synthesis for Online Optimisation of Batch Reactors with Multiple Manipulated Inputs. Comput. Chem. Eng. 1998, 22 10, 1429. (2) Knopf, F. C.; Okos, M. R.; Reklaitis, G. V. Optimal Design of Batch/Semicontinuous Processes. Ind. Eng. Chem. Process Des. Dev. 1982, 21, 79. (3) Yeh, N. C.; Reklaitis, G. V. Synthesis of Batch/Semicontinuous Processes. Comput. Chem. Eng. 1987, 11, 639. (4) Ravemark, D. E.; Rippin, D. W. T. Optimal Design of Multiproduct Batch Plants. Comput. Chem. Eng. 1998, 22, 177. (5) Montagna, J. M.; Vecchietti, A. R.; Iribarren, O. A.; Pinto, J. M.; Asenjo J. A. Optimal Design of Protein Production Plants with Time and Size Factor Process Models. Biotechnol. Prog. 2000, 16, 228. (6) Salomone, H. E.; Iribarren, O. A. Posynomial Modeling of Batch Plants. Comput. Chem. Eng. 1992, 16, 173. (7) Pinto, J. M.; Montagna, J. M.; Vecchietti, A. R.; Iribarren, O. A.; Asenjo, J. A. Process Performance Models in the Optimiza-
tion of Multiproduct Protein Production Plants. Biotechnol. Bioeng. 2001, 74, 451. (8) Barrera, M. D.; Evans, L. B. Optimal Design and Operation of Batch Processes. Chem. Eng. Commun. 1989, 82, 45. (9) Salomone, H. E.; Montagna, J. M.; Iribarren, O. A. Dynamic Simulations in the Design of Batch Processes. Comput. Chem. Eng. 1994, 18, 191. (10) Bathia, T. K.; Biegler, L. T. Dynamic Optimization in the Design and Scheduling of Multiproduct Batch Plants. Ind. Eng. Chem. Res. 1996, 35, 2234. (11) Uesbeck, F.; Samsatli, N. J.; Papageorgiou, L. G.; Shah, N. Robust Optimal Fermentation Operating Policies. Comput. Chem. Eng. 1998, 22, S167. (12) Monod, J. The Growth of Bacterial Cultures. Annu. Rev. Microbiol. 1949, 3, 371. (13) De la Cruz Soriano, R.; Gonza´lez Sua´rez, E.; Pedraza Ga´rciga, J. Produccio´ n de Acetal a partir de Bioetanol, 1st ed.; CYTED, Ciencia y Tecnologı´a para el Desarrollo: Buenos Aires, Argentina, 2003. (14) Atkinson, K. E. An Introduction to Numerical Analysis, 2nd ed.; John Wiley & Sons: New York, 1989. (15) Brooke, A.; Kendrick, D.; Meeraus, A.; Raman, R. GAMS, A User’s Guide; Scientific Press: Redwood City, CA, 1998. (16) Oh, M.; Pantelides, C. C. A Modelling and Simulation Language for Combined Lumped and Distributed Parameter Systems. Comput. Chem. Eng. 1996, 20, 611.. (17) Barton, P. I.; Pantelides, C. C. Modelling of combined discrete/continuous processes. AIChE J. 1994, 40, 966. (18) Asenjo, J. A.; Montagna, J. M.; Vecchietti, A. R.; Iribarren, O. A.; Pinto, J. M. Strategies for the Simultaneous Optimisation of the Structure and the Process Variables of a Protein Production Plant. Comput. Chem. Eng. 2000, 24, 2277. (19) Petrides, D.; Sapidou, E.; Calandranis, J. Computer-Aided Process Analysis and Economic Evaluation for Biosynthetic Human Insulin ProductionsA Case Study. Biotechnol. Bioeng. 1995, 48, 529.
Received for review July 2, 2003 Revised manuscript received May 6, 2004 Accepted May 17, 2004 IE030549H
