Dynamic Optimization of a Continuous Polymer Reactor Using a

The proposed modified differential evolution (MDE) algorithm is different from differential evolution algorithms in the sense that MDE employs a local...
0 downloads 0 Views 89KB Size
Download PDF
Ind. Eng. Chem. Res. 1999, 38, 4825-4831

4825

Dynamic Optimization of a Continuous Polymer Reactor Using a Modified Differential Evolution Algorithm Moo Ho Lee, Chonghun Han,* and Kun Soo Chang Automation Research Center and Department of Chemical Engineering, Pohang University of Science and Technology, Pohang 790-784, Korea

Dynamic optimization of a continuous polymer reactor aims to decide optimal trajectories of control input variables so that the transition time, required to reach the desired steady state from the initial state during startup or grade-change operation, is minimized. The problem is challenging because of its highly nonlinear dynamics and multimodal properties. The proposed modified differential evolution (MDE) algorithm is different from differential evolution algorithms in the sense that MDE employs a local search to enhance the computational efficiency and modified heuristic constraints to systematically reduce the size of the search space. The algorithm is illustrated by several case studies using the dynamic optimization problem of a continuous methyl methacrylate-vinyl acetate copolymerization reactor. The case studies have shown that the proposed algorithm offers faster speed, flexible implementation, and higher robustness to find the global optimum than differential evolution algorithms. 1. Introduction Batch process operations are popular ways of producing high value-added products such as pharmaceutical and specialty chemical products. The goal of batch process operators is to decide the optimal trajectories of input variables so that the residence time can be minimized and the qualities of the product are within the product specifications. The dynamic optimization (optimal control) problem has been formulated and solved to find optimal control input trajectories in a batch or a semibatch reactor since the 1950s.1-5 Dynamic optimization problems of continuous reactors during the startup or grade-change operation periods are similar to the dynamic optimization of a batch reactor. The goal of the startup operation for a continuous reactor is to operate the process so that the desired steady state is reached from the initial state within the shortest time. On the other hand, the goal of the grade-change operation is to change within the shortest time the state of the process from one condition to another where product specifications for new grades of the product are satisfied. The transition time, which it takes for the process to reach the goal state from the initial state, is economically very important because the products produced during the transition time are offspecification and sent to waste treatment facilities.6 The calculus of variation-based optimal control methods has been extensively applied to the dynamic optimization of a batch or semibatch.3,7-9 However, the calculus of variation-based optimal control method is efficient only for simple optimal control problems because of its extensive computation cost. Iterative dynamic programming (IDP), which is based on the principle of optimality in the dynamic programming, has been proposed by Luus10,11 and applied to large-size problems successfully. The third approach is the parametrization method, where the state or control input trajectory is discretized * Author to whom correspondence should be addressed. E-mail: [email protected]. Fax: (+82) 562-279-3499. Tel: (+82) 562-279-2279.

and approximated as several functions. Complete parametrization methods12,13 take into account state constraints. However, they are more difficult to implement and solve than the control parametrization methods.14-17 In general, because chemical reaction systems are nonlinear, multimodal, and discontinuous, both complete and control parametrization methods have difficulties in the convergence to the global optimum.18 For example, the use of control parametrization methods with deterministic optimization algorithms based on gradient information, such as successive quadratic programming (SQP), frequently gives a local optimum as a result. Thus, stochastic optimization methods such as genetic algorithm (GA)19-22 and evolution algorithm (EA) strategies23 have been studied as promising candidates for the optimization of discontinuous and multimodal chemical systems. Although EAs are easy to implement and efficient in finding the global optimum, they usually need more computation time than deterministic optimization methods because of the use of probabilistic transition rules. Lee et al.24 have proposed a two-level hierarchical dynamic optimization methodology for continuous processes by combining heuristic constraints with genetic algorithms to reduce the computation time. The differential evolution (DE), which was first proposed by Storn and Price,25 is very simple to understand and implement. Furthermore, it needs only a few control parameters. Wang and Chiou26 have used DE to generate initial starting points for SQP in the optimal control and the optimal time location problems of differential-algebraic systems. To apply DE to a problem which needs a significant amount of computation time in the evaluation of the objective function such as dynamic optimization, the convergence speed of DE has to be improved sufficiently to dramatically reduce the number of objective function evaluations. Chiou and Wang27 proposed a hybrid differential evolution algorithm and applied it to the parameter estimation problem. The algorithm uses a descent search to improve the performance of DE.

10.1021/ie980373x CCC: $18.00 © 1999 American Chemical Society Published on Web 11/06/1999

4826

Ind. Eng. Chem. Res., Vol. 38, No. 12, 1999

In this paper, we propose a modified differential evolution (MDE) algorithm for the dynamic optimization problem during the startup or grade-change period of continuous processes. MDE has two characteristics that are different from those of DE. First, it uses a local search based on the difference vector information that enhances the computation speed without loss of the robustness of DE. Second, the algorithm employs modified heuristic constraints called converging range constraints to systematically reduce the size of the search space.24 2. Dynamic Optimization Using MDE 2.1. Dynamic Optimization Problem. The dynamic optimization problem for a continuous process is to find the control input vector u(t) over t ∈ [t0, tf] that minimizes the performance index J(u(t),x(t)). This problem can be formulated as:

Min J(u(t),x(t)) ) u(t)

∫tt g(x(t),u(t),t) dt f

0

(1)

for control input when converging range constraints are applied and are obtained as follows: R ui,j ) uss i -

(u j i - ui) j Ri 2

(9)

R u j i,j ) uss i +

(u j i - ui) j Ri 2

(10)

R R ui,j g ui and u j i,j eu ji

(11)

Rji is the jth power of converging parameter Ri for the ith control input. According to the above equations, Rji R R goes close to zero and ui,j and u j i,j approach uss i , as j increases. and The upper bound for the transition time tmax f steady-state control inputs uss i to meet the specified product properties after startup or grade change are found at the upper level of a two-level hierarchical dynamic optimization with converging range constraints for the continuous processes. The upper level optimization problem can be formulated as follows:

Min J(u(t),x(t)) )

subject to

max uss,tf

dx/dt ) f(x(t),u(t),t), x(t0) ) x0

(2)

ui e ui(t) e u j i, i ) 1, ..., m

(3)

u(t),R,S

∫0t

max f

g(x(t),u(t),t) dt

max f

0

g(x(t),uss,t) dt +

∫tt

∞ f max f

g(x(t),uss,t) dt (12)

subject to

where x(t) and u(t) denote the state and the control input, respectively. ui and u j i are the minimum and maximum bounds for the control variable ui. t0 is the starting time for grade change or startup, and tf is a sufficiently large time for the process to reach the steady state. For most of continuous processes, g(x(t),u(t),t) represents the deviation from the desired product properties. Thus, if we choose tf as a sufficiently large value and minimize the objective function, the desired steady state will be achieved and maintained.15 The transition time [t0, tf] is divided into n time intervals to parametrize the control input values at each time interval. Each ui,j, i ) 1, ..., m, j ) 1, ..., n, is assumed to have a constant value during the [tj-1, tj] time interval. The two-level hierarchical dynamic optimization for continuous processes24 has been introduced to reduce the size of the search space and enhance the search speed without loss of search robustness. The dynamic optimization problem can be reformulated based on the two-level hierarchical dynamic optimization as follows:

Min J(u(t),x(t)) )

∫tt

(4)

subject to

dx/dt ) f(x(t),u(t),t), x(t0) ) x0

(5)

R R e ui,j e u j i,j , i ) 1, ..., m, j ) 1, ..., n ui,j

(6)

0eRe1

(7)

Si e Si,j e Si, i ) 1, ..., m, j ) 1, ..., n

(8)

R R ui,j and u j i,j are the minimum and the maximum bounds

dx/dt ) f(x(t),u(t),t), x(t0) ) x0

∫tt

∞ f max f

g(x(t),uss,t) dt < 

(13) (14)

ss j ss uss i e ui e u i , i ) 1, ..., m

(15)

tmax e tmax e ht max f f f

(16)

If the desired steady-state control input values are applied as the constant inputs to the process during the transition time, the process will reach the desired steady . This state after the maximum transition time tmax f is used as the upper maximum transition time tmax f bound for the transition time at the low level where both optimal trajectories of control input variables and the optimal transition time will be found. Equation 14 is the condition used to check whether the process has reached the steady state. t∞f is sufficiently large time to reach the steady state. To decide the optimal switching time of control inputs, the time interval size parameter Si,j ) ti,j - ti,j-1 has been introduced for each interval [ti,j - ti,j-1]. Si,j also has the minimum and the maximum bounds like the control input.15,24,26,28 2.2. MDE Algorithm. MDE has been developed based on DE for the dynamic optimization of continuous processes. MDE has the following two features: (a) search performance improvement from a local search and (b) systematic reduction of search space from the introduction of heuristic constraints. The better search performance of DE than conventional genetic algorithms comes from the mutation operation that uses a difference vector. A difference vector is composed of two randomly chosen individuals and plays the role of search direction in DE like a gradient vector in the deterministic optimization. The optimal step size that determines

Ind. Eng. Chem. Res., Vol. 38, No. 12, 1999 4827

the step length in the direction of the difference vector should be decided to get the best search performance. However, DE uses a constant value for the step size called the weighting factor F during the entire search process. In MDE, we evaluate the influence of the difference vector on the objective function and utilize the information to find a better step size through a simple local search. If the newly generated individual has a better objective function value than the current one, then the weighting factor F is increased and the new individual will be generated using the new increased weighting factor F and the same difference vector. This process continues until there is no more improvement in the objective function value during the local search. The converging range constraints24 are incorporated into MDE more easily and efficiently. In MDE, the converging parameter R, control input variables, and time interval sizes are optimized simultaneously to improve search speed. All of the steps of the MDE algorithm are as follows: First, we define the optimization parameter vector z that is composed of control vector u, time interval size vector S, and converging parameter R. Step 1. Choose a population size NP, a base weighting factor Fb, a crossover factor CR, and a step factor ξ. The weighting factor is F ∈ [0, 1] in DE. As the local search of MDE reduces search robustness, the base weighting factor Fb is given a smaller value than the weighting factor value used in DE to enhance the robustness of the search in MDE. Step 2. Set F :) Fb. Step 3. Generate the initial population by generating random numbers βk.

zk ) z + βk(zj - z), k ) 1, ..., NP

(17)

R where zk ) [ui,j Si,j Ri, i ) 1, ..., m, j ) 1, ..., n]k, z ) [ui,j R Si,j Ri, i ) 1, ..., m, j ) 1, ..., n], and zj ) [u j i,j S h i,j R j i, i ) 1, ..., m, j ) 1, ..., n] are the kth individual, minimum, and maximum bound vectors for the optimization parameter vector z, respectively. βk represents uniformly distributed random numbers. Step 4. Perform the mutation operation using two randomly selected individuals and the best-so-far individual, zbG at the Gth generation.

G G G zˆ G+1 ) zkG + F(zG k b - zk ) + F(zr1 - zr2),

k ) 1, ..., NP and k * r1, r2 (18) where r1 and r2 are randomly determined to select two individuals and F ∈ [0, 1] is a weighting factor. We will obtain faster convergence and less robust results as F has a larger value.26 Several strategies have been proposed for the mutation operation.25,29 For example, a difference vector may be composed of either two randomly chosen individuals or a pair of the best-sofar individuals along with the two randomly chosen ones. A new individual can be generated using one difference vector or the combination of two difference vectors. As these strategies have shown somewhat different performances depending upon the problems, the above mutation strategy has been chosen for our problem through simulation runs.

Step 5. Perform a crossover operation to increase the offspring diversity as G+1 zk,l )

{

G+1 zk,l CR e δk,l G+1 zˆ k,l CR > δk,l k ) 1, ..., NP and l ) 1, ..., L (19)

where l is the index for the optimization parameter of zk, δk,l is the uniformly distributed random number, and the crossover factor CR ∈ [0, 1] determines the diversity of generated offsprings.29 The larger the value CR has, G+1 the higher the probability of being an offspring zk,l G+1 the newly generated individual zˆ k,l will have.26 Step 6. Check whether newly generated individuals are within their bounds. If the value of an individual is out of bounds, replace it with its minimum or its maximum value. G+1 ) zk,l

{

G+1 zl, zk,l < zl G+1 zjl, zk,l > zjl k ) 1, ..., NP and l ) 1, ..., L (20)

Step 7. Check whether a local search is needed. ), k ) 1, (i) Evaluate the objective function J(zG+1 k ..., NP, for each individual. ) < J(zG (ii) If J(zG+1 k k ), then go to step 8. Else go to step 9. Step 8. Perform the local search by adjusting the weighting factor F. (i) F ) ξF where ξ is a step factor that controls the next search step size. Because the base weighting factor value is smaller than the weighting factor value used in DE, the weighting factor value is increased gradually by increasing the step factor. The range for the step factor, ξ ∈ [1.3, 1.7], has been determined from simulation runs. :) zG+1 . (ii) zG-new k k through eqs (iii) Generate a new individual zG+1 k 18-20 and the new F. At this step, we use the same G G G values for zG k , zb , zr1, and zr2 as the ones at the previous generation in eq 18 in order that the new individual has the same search direction as the previous one. (iv) If zG+1 < zG-new , then go back to step 8(i). k k :) zG-new and go to step 9. Otherwise, zG+1 k k Step 9. Select the individuals of the next generation . and the best individual zG+1 b G G+1 ) < J(z ), then z :) zG+1 . If J(zG+1 k k k k G+1 G Otherwise, zk :) zk . Step 10. Set F :) Fb. Step 11. If the maximum generation number or the desired objective function value is obtained, then stop. Otherwise, go back to step 4. 3. Dynamic Optimization of a Methyl Methacrylate-Vinyl Acetate Copolymerization Reactor The performance of MDE for dynamic optimization of continuous processes is illustrated by its application to the dynamic optimization problem of a continuous methyl methacrylate-vinyl acetate (MMA-VA) copolymerization reactor.30 The free radical solution MMAVA copolymerization model has been chosen due to its complicated and highly nonlinear characteristics.

4828

Ind. Eng. Chem. Res., Vol. 38, No. 12, 1999

The objective of this dynamic optimization problem is to determine the optimal trajectories of control input values that minimize the transition time. After the transition time, the process should reach the desired steady state and the product qualities should satisfy product specifications. Controlled variables are the weight-average molecular weight (Mw) and the monomer composition of VA in dead copolymer (YVA). The jacket temperature (K), u1, and the flow rate of the monomer VA (kg/h), u2, have been chosen as control input variables. The control objective is

Min J ) u1,u2

∫tt ($1|Mw(t) - Mwd| + f

0

$2|YVA(t) - YVAd|) dt 298.15 e u1 e 373.15 and 0 e u2 e 200 (21) where $1 ) 1.0 × 10-7 and $2 ) 1.0 × 10-2. The desired molecular weight (Mwd) is 50 000 and the desired monomer composition of VA in dead copolymer (YVAd) is 0.2. To show the search performance of MDE, we compare MDE with DE26 in section 3.1. The optimal switch time for control inputs has been decided in section 3.2. 3.1. Comparisons of MDE with DE. At the upper level of two-level hierarchical dynamic optimization for has been chosen as 150 000 continuous processes, tmax f ss and u have been determined as 326.2 (K) s and uss 1 2 and 52.4 (kg/h), respectively. The minimum and maxiare 50 000 s and 200 000 s, mum bounds for tmax f respectively, and the tolerance for steady state, , is 3.0. ] into 15 time At the low level, we divide [t0, tmax f intervals and define control input parameters, ui,j, i ) 1, 2 and j ) 1, ..., 15, for 15 time intervals. The control parameters for DE and MDE have been chosen from many simulation runs based on the values recommended by Storn,29 so that DE and MDE could show their best performances. We have chosen NP ) 50, CR ) 0.8, and F ) 0.7 for DE and NP ) 50, CR ) 0.8, Fb ) 0.5, and step factor ξ ) 1.5 for MDE. The value of the base weighting factor Fb in MDE is usually determined as the smaller value than that of the weighting factor F in DE to increase the robustness of MDE. The number of local searches in MDE was observed to be about 300500 during 50 000 objective function evaluations. Figure 1 shows the results of DE proposed by Wang and Chiou,26 DE with only a local search, DE with only converging range constraints, and MDE. MDE gives faster convergence and a better objective function value than DE. Wang and Chiou26 used DE to generate the best initial starting points for SQP and then applied the SQP solver to obtain a more accurate control input trajectory. They also used a modified collocation method to solve the differential algebraic equations (DAEs). As our copolymerization reactor has complex and very nonlinear ordinary differential equations (ODEs), the Adams-Gear method, which is suitable for stiff systems such as polymerization reactors, has been used to solve them.31 SQP has been applied to our problems with the initial starting points which were obtained from DE with 1000 and 5000 objective function evaluations. However, it failed to find the better solutions for the two cases. Thus, only DE has been applied to our problem without SQP. It is shown in the results of DE and DE with converging range constraints that the converging range

Figure 1. Comparisons of objective function values with the number of objective function evaluations among DE (Wang and Chiou), DE with only local search, DE with only converging range constraints, and MDE.

Figure 2. Converging parameter R of the best individual according to the number of objective function evaluations during the optimization in MDE.

constraints enable the faster convergence of optimization. When DE and DE with local search are compared, it is easily shown how local search in MDE is efficient in the improvement of search speed without the loss of search robustness. The objective function values of DE and MDE are 65.2 and 16.2, respectively, at 5000 objective function evaluations. For the comparisons of MDE with other methods, the objective function values of GAs, IDP, and SQP are 67.1, 50.7, and 89.6, respectively, at 5000 objective function evaluations.24 Figure 2 shows the converging parameter R of the best individual according to objective function evaluations during the optimization in MDE. As shown in Figure 2, a smaller value of R enables the rapid convergence of MDE at the beginning of optimization and a larger value of R allows MDE to investigate a wider search area at the latter part of optimization. When R and control input trajectories are simultaneously optimized, a better performance in terms of search speed and objective function value is obtained. Figures 3-6 show

Ind. Eng. Chem. Res., Vol. 38, No. 12, 1999 4829

Figure 3. Molecular weight trajectories for DE (Wang and Chiou) and MDE after 5000 objective function evaluations.

Figure 5. Optimal jacket temperature policy (u1) for DE (Wang and Chiou) and MDE after 5000 objective function evaluations.

Figure 4. VA composition trajectories for DE (Wang and Chiou) and MDE after 5000 objective function evaluations.

Figure 6. Optimal VA feed flow rate policy (u2) for DE (Wang and Chiou) and MDE after 5000 objective function evaluations.

the obtained output and control trajectories after 5000 objective function evaluations in DE and MDE. 3.2. Optimal Time Interval Size. To investigate the effect of time interval size on the dynamic optimization of a continuous copolymerization reactor, the optimal time interval sizes for all of the discretized time intervals have been determined using MDE. In equal time intervals, the control inputs u1 and u2 have an equal time interval of 10 000 s for all time intervals and all control variables. Therefore, no time size parameter Si,j is needed. In variable time intervals, u1 and u2 have different time interval size parameters, S1,j and S2,j, respectively, at jth time interval. Thus, the number of time interval size parameters becomes 30 and time interval size parameter Si,j has a value between 100 and 20 000 s. The results of dynamic optimization are shown in Table 1. The variable time size interval approach shows the better results in terms of the objective function value than the equal time size interval method because of a more accurate switching time. Figures 7 and 8 show the optimal control input trajectories for equal time intervals and variable time

Table 1. Comparisons of Equal Time Intervals and Variable Time Intervals no. of objective function evaluations 5000 50 000

objective function value no. of local searches CPU time (s) objective function value no. of local searches CPU time (s)

equal time intervals

variable time intervals

16.24 284 378 15.01 457 4098

7.08 328 417 4.78 532 4297

intervals using MDE. The variable time size interval approaches have a shorter or longer switching time for each time interval to take into account more accurate switching time. 4. Conclusions A MDE algorithm for the dynamic optimization of continuous processes has been proposed. MDE has been applied to the dynamic optimization of a continuous MMA-VA copolymerization reactor to decide the opti-

4830

Ind. Eng. Chem. Res., Vol. 38, No. 12, 1999

Figure 7. Comparisons of obtained optimal jacket temperature policy (u1) between equal time intervals and variable time intervals after 5000 objective function evaluations in MDE.

Fb ) base weighting factor of MDE J ) objective function Mw ) molecular weight Mwd ) desired molecular weight NP ) number of population in DE, MDE r1, r2 ) random number S ) interval size parameter t ) time t0 ) initial time tf ) final time t∞ ) sufficiently large time after transition time tmax ) upper bound for transition time f , ht max ) minimum and maximum bound of the upper tmax f f bound for transition time u ) control input vector u, u j ) minimum and maximum bound of the control input vector j R ) minimum and maximum bound of the control uR, u input vector in converging range constraints u1 ) jacket temperature (K) u2 ) feed flow rate of monomer VA (kg/h) uss ) control input for the steady state x ) state vector YVA ) monomer composition of VA in the dead copolymer YVAd ) desired monomer composition of VA in the dead copolymer z ) optimization parameter vector z, zj ) minimum and maximum bound of the optimization parameter vector Subscripts and Superscripts G ) index for generation in MDE i ) index for the control variable j ) index for the time interval k ) index for the individual in MDE l ) index for the optimization parameter Greek Symbols

Figure 8. Comparisons of obtained optimal VA feed flow rate policy (u2) between equal time intervals and variable time intervals after 5000 objective function evaluations in MDE.

mal trajectories for input variables that can minimize the transition time during startup or grade-change operation. Case studies have shown that the local search and the heuristic constraints introduced for the reduction of search space size have allowed MDE to have the faster speed and the better performance in terms of objective function values than DE without loss of the robustness of evolution algorithms. Acknowledgment The authors thank R. Luus, J. P. Chiou, and F. S. Wang for their valuable comments and support and the Korea Science and Engineering Foundation (KOSEF) for financial support through Automation Research Center at POSTECH. The authors also thank the school of Environmental Engineering at POSTECH for partial support. Nomenclature CR ) crossover factor of DE and MDE F ) weighting factor of DE and MDE

R ) converging parameter β ) uniformly distributed random number δ ) uniformly distributed random number  ) tolerance for the steady state ω1, ω2 ) weighting factors of the objective function ξ ) step factor

Literature Cited (1) Denbigh, K. G. Optimum Temperature Sequences in Reactors. Chem. Eng. Sci. 1958, 8, 125. (2) Edgar, T. F.; Lapidus, L. The Computation of Optimal Singular Bang-Bang Control II: Nonlinear Systems. AIChE J. 1972, 18, 780. (3) Wu, G. Z. A.; Denton, L. A.; Laurence, R. L. Batch Polymerization of Styrene-Optimal Temperature Histories. Polym. Eng. Sci. 1982, 22, 1. (4) San, K. Y.; Stephanopoulos, G. Optimization of a Fed-Batch Penicillin Fermentation Processes. Biotechnol. Bioeng. 1989, 34, 72. (5) Lee, J.; Ramirez, W. F. On-Line Optimal Control of Induced Foreign Protein Production by Recombinant Bacteria in Fed-Batch Reactors. Chem. Eng. Sci. 1996, 51, 521. (6) Sirohi, A.; Choi, K. Y. On-Line Parameter Estimation in a Continuous Polymerization Process. Ind. Eng. Chem. Res. 1996, 35, 1332. (7) Ponnuswamy, S. R.; Shah, S. L.; Kiparissides, C. A. Computer Optimal Control of Batch Polymerization Reactor. Ind. Eng. Chem. Res. 1987, 26, 2229. (8) Takamatsu, T.; Shioya, S.; Okada, Y. Molecular Weight Distribution Control in a Batch Polymerization Reactor. Ind. Eng. Chem. Res. 1988, 27, 93.

Ind. Eng. Chem. Res., Vol. 38, No. 12, 1999 4831 (9) Chang, J. S.; Lai, J. Computation of Optimal Temperature Policy for Molecular Weight Control in a Batch Polymerization Reactor. Ind. Eng. Chem. Res. 1992, 31, 861. (10) Luus, R. Optimal Control by Dynamic Programming using Accessible Grid Points and Region Reduction. Hung. J. Ind. Chem. 1989, 17, 523. (11) Luus, R. Application of Dynamic Programming to HighDimensional Nonlinear Optimal Control Problems. Int. J. Control 1990, 52, 239. (12) Biegler, L. T. Solution of dynamic optimization problems by successive quadratic programming and orthogonal collocation. Comput. Chem. Eng. 1984, 8, 243. (13) Cuthrell, J. E.; Biegler, L. T. Simultaneous optimization and solution methods for batch reactor control profiles. Comput. Chem. Eng. 1989, 13, 49. (14) Goh, C. J.; Teo, K. L. Control Parameterization: a Unified Approach to Optimal Control Problems with General Constraints. Automatica 1988, 24, 3. (15) McAuley, K. B.; MacGregor, J. F. Optimal Grade Transitions in a Gas Phase Polymer Reactor. AIChE J. 1992, 38, 1564. (16) Vassiliadis, V. S.; Sargent, R. W.; Pantelides, C. C. Solution of a Class of Multistage Dynamic Optimization Problems. 1. Problems without Path Constraints. Ind. Eng. Chem. Res. 1994, 33, 2111. (17) Barton, P. I.; Allgor, R. J.; Feehery, W. F.; Galan, S. Dynamic Optimization in a Discontinuous World. Ind. Eng. Chem. Res. 1998, 37, 966. (18) Luus, R. Optimal Control of Batch Reactors by Iterative Dynamic Programming. J. Process Control 1994, 4, 218. (19) Goldberg, D. E. Genetic Algorithms in Search, Optimization and Machine Learning; Addison-Wesley: New York, 1989. (20) Michalewicz, Z.; Janikow, C. Z.; Krawczyk, J. B. A Modified Genetic Algorithm for Optimal Control Problems. Comput. Math. Appl. 1992, 23, 83. (21) Seywald, H.; Kumar, R. R.; Deshpande, S. M. Genetic Algorithm Approach for Optimal Control Problems with Linearly Appearing Controls. J. Guidance, Control, Dynamics 1995, 18, 177.

(22) Park, S.; Yoon, E. S. Study on Chemical Process Synthesis and Optimization Strategies using Genetic Algorithms. Proc. ’96 KIChE Fall Meeting 1996, 2, 1463. (23) Michalewicz, Z. Genetic Algorithms + Data Structures ) Evolution Programs; Springer-Verlag: New York, 1996. (24) Lee, M. H.; Han, C.; Chang, K. S. Hierarchical TimeOptimal Control of a Continuous Copolymerization Reactor during Start-up or Grade Change Operation using Genetic Algorithms. Comput. Chem. Eng. 1997, 21, S1037. (25) Storn, R.; Price, K. Minimizing the Real Functions of the ICEC’96 Contest by Differential Evolution. Proceedings of the IEEE Conference on Evolutionary Computation, Nagoya, Japan, 1996; p 842. (26) Wang, F. S.; Chiou, J. P. Optimal Control and Optimal Time Location Problems of Differential-Algebraic Systems by Differential Evolution. Ind. Eng. Chem. Res. 1997, 36, 5348. (27) Chiou, J. P.; Wang, F. S. Hybrid Differential Evolution for Parameter Estimation of a Batch Bioprocess. Proceedings of the IEEE Symposium on Control Theory and Applications, Singapore, 1997; p 171. (28) Bojkov, B.; Luus, R. Time Optimal Control of Nonlinear Systems with Unspecified Final Times. Chem. Eng. Sci. 1996, 51, 905. (29) Storn, R. On the Usage of Differential Evolution for Function Optimization. NAFIPS, Berkeley, 1996; p 519. (30) Congalidis, J. P.; Richards, J. R.; Ray, W. H. Feedforward and Feedback Control of a Solution Copolymerization Reactor. AIChE J. 1989, 35, 891. (31) Gear, C. W. Numerical Initial Value Problems in Ordinary Differential Equations; Prentice-Hall: Englewood Cliffs, NJ, 1971.

Received for review June 9, 1998 Revised manuscript received August 12, 1999 Accepted September 2, 1999 IE980373X