Trajectory Tracking of an Optimizing Path in a Batch Reactor

Jun 1, 1996 - trajectory. Yield optimization is assured for a real plant by tracking the model-dependent optimized trajectory through the proposed mod...
0 downloads 0 Views 294KB Size
Ind. Eng. Chem. Res. 1996, 35, 2247-2260

2247

Trajectory Tracking of an Optimizing Path in a Batch Reactor: Experimental Study Jyh-Shyong Chang,* Jen-Sheng Hsu, and Yao-Tsung Sung Department of Chemical Engineering, Tatung Institute of Technology, Taipei, Taiwan, ROC

We propose an integrated method for optimization and control of batch reactors. Based on the desired performance index, the modified two-step method is applied to optimize an operating trajectory. Yield optimization is assured for a real plant by tracking the model-dependent optimized trajectory through the proposed modified globally linearizing control (MGLC) structure. Experimental results reveal that the proposed MGLC structure can be applied in tracking an operating trajectory determined on-line or off-line. 1. Introduction There is considerable interest in the control of batch chemical reactors (Jutan and Uppal, 1984; Liptak, 1986; Jutan and Rodriguez, 1987; Kravaris and Chung, 1987; Kravaris et al., 1989; Lewin and Lavie, 1990; Chang and Huang, 1994; Chang and Chung, 1995). From an engineering point of view, the control of batch chemical reactors is important for optimization of the yield (Bhat et al., 1991) and the quality of a product (Ponnuswamy and Shah, 1987; Soroush and Kravaris, 1992; Ellis et al., 1994). Both problems can be formulated as tracking problems: given

x3 ) f(x,u) y ) h(x)

(1)

The objective is to find a control law so that y(t) tracks a given trajectory y*(t) (Kravaris et al., 1989). In this matter, whether one can operate a batch reactor efficiently hinges on proper choice of the set-point policy and the control structure for trajectory tracking. Many simulation studies have proposed control strategies in dealing with the trajectory-tracking problem of batch reactors (Jutan and Uppal, 1984; Liptak, 1986; Jutan and Rodriguez, 1987; Kravaris and Chung, 1987; Kravaris et al., 1989; Lewin and Lavie, 1990). Our previous works (Chang and Huang, 1994; Chang and Hseih, 1995; Chang and Chung, 1995) were also devoted to dealing with the trajectory-tracking problem by proposing the modified globally linearizing control (MGLC) structure. By simulating different exothermic batch and semibatch reactors in which a simple chemical reaction or polymerization reactions were carried out, we had tested the MGLC structure over a wide range of system parameters. The behavior of the proposed MGLC structure was proven to be predictable and reliable with tuning parameters based on the proposed tuning method if the manipulated inputs are not constrained. To support the validity of the MGLC structure in the real world, this study aims to set up a batch reactor system and test the proposed MGLC structure on it. The batch reactor consists of a conventional vessel equipped with a heating jacket and cooling * Author to whom correspondence should be addressed. E-mail: [email protected]. FAX: 886-2-5861939.

S0888-5885(95)00534-3 CCC: $12.00

coils. We designed the experimental apparatus in this way because the heating and cooling systems designed can be constructed and scaled up much easier than other arrangements (Soroush and Kravaris, 1992). Moreover, the set-point policy for batch reactors has received much less attention in connection with the geometric nonlinear control. Palanki et al. (1994) synthesized the optimal state feedback laws for the yield optimization problem. The state inequality constraints and free terminal time were considered in his works. Based on the maximum principle and the differential geometric approach, the end-point optimization of the free time problem was treated successfully (Palanki et al., 1994). However, implementation of the proposed optimal state feedback laws adopting full dynamics of the system in a real process involves tackling the problems of model uncertainties and parameter variations. If there exists any model/plant mismatch, the control target will be missed. From the engineering viewpoint, it may be much easier to obtain accurate kinetic data alone for a chemical reaction. Instead, identifications of the performance of a batch reactor system such as heat transfer, agitation, and metal dynamics are difficult. Therefore, there would be more incentive to adopt the reliable kinetic model alone than the integrated dynamic model of the system in developing a set-point policy for batch reactors. Another comment on the work of Palanki et al. (1994) is that a combination of the geometric approach (Isidori, 1989) and Pontryagin’s principle (Sage and White, 1977) may not be easily applicable especially when the analytical differentiation of the system equations is difficult to manipulate. On the other hand, because the computation algorithm of the modified two-step method is based on integration of the system equations only rather than on differentiation manipulation, adoption of the modified two-step method for off-line or on-line determination of an optimal operating trajectory is attractive. For online computation, this method uses the states provided by the state observer usually built in the MGLC structure (Chang and Chung, 1995). In this way, setpoint optimization can be linked with the geometric nonlinear control effectively. The modified two-step method developed in our previous researches emphasized tackling the product quality optimization problem (Takamatsu et al., 1988; Chang and Lai, 1992; Chang and Chung, 1995); therefore, in © 1996 American Chemical Society

2248

Ind. Eng. Chem. Res., Vol. 35, No. 7, 1996

this study, we will examine the applicability of the modified two-step method to the yield optimization problem occurring in a batch reactor. This article first proposes a set-point policy for yield optimization by extending the modified two-step method (Chang and Lai, 1992; Chang and Chung, 1995). Then a brief review of the MGLC structure (Chang and Huang, 1994; Chang and Hseih, 1995; Chang and Chung, 1995) is presented. Application of the proposed optimization and control method will be implemented in tracking an optimizing temperature trajectory of a batch reactor. Synthesis of hexyl monoester of maleic acid is the reaction studied in this work (Hugo, 1981; Westerterp et al., 1984). 2. Calculation of an Optimizing Path by the Modified Two-Step Method For the batch reactor control system, the governing equations excluding the one that contains the actual manipulated input u in eq 2 are usually formulated on the kinetic model of the reaction under consideration. Therefore, the general dynamics for batch reactors can usually be expressed as

dx1 ) f1(x) dt

dx2 ) f2(x1,x2,...,xn-1,u′) dt l dxn-1 ) fn-1(x1,x2,...,xn-1,u′) dt

(3)

To simplify the discussion, a small system with n ) 3 is now considered. A similar conclusion can be deduced for the system with n > 3. For the yield optimization problem, the performance indices are usually chosen as a function of the terminal values of x1, x2, and u′ at tf. Then whether the desired performance index can be achieved hinges on how the following internal state trajectory is tracked:

dx2 f2(x1,x2,u′) ) dx1 f1(x1,x2,u′) ) ξ(x1)

(4)

where x1 ∈ [x1(t0), x1(tf)]. Any three parameter solution forms of ξ(x1) can be adopted (Chang and Lai, 1992). In this study, we choose a second order polynomial.

dx2 ) f2(x) dt

ξ(x1) ) a0′ + a1′x1 + a2′x12

(5)

The parameters a0′, a1′, and a2′ are to be determined by solution of three simultaneous algebraic equations, which are formulated by the chosen solution form and the following three conditions proposed in this study:

l dxn-1 ) fn-1(x) dt dxn ) fn(x,u) dt

dx1 ) f1(x1,x2,...,xn-1,u′) dt

(2)

Among the subsystems of eq 2 excluding the last differential equation, we can choose one state variable as the optimizing variable u′. The chosen optimizing variable shall have a strong and direct effect on the performance indices for operating a batch reactor and can be easily manipulated. Nonisothermal temperature trajectories are a common means for optimizing the throughput of many batch reactors (Soroush and Kravaris, 1993). Here, the performance indices for the yield optimization problem are usually chosen to be a function of terminal states. In essence, a suitable choice of terminal states can define the optimum for operating batch reactors, such as least end time, operational safety, and desired yield. Therefore, application of the proposed modified two-step method for batch reactors is to find an optimizing variable, u′(t) ) xn(t), to drive the states of the reactor along a suitable internal state trajectory determined by the boundary conditions. When the whole internal state trajectory can be tracked properly, the desired terminal states can be achieved. Because the modified two-step method adopts the sufficient conditions for optimality, the optimum can be guaranteed if the calculated optimizing variable, u′(t), is exact and physically possible along the trajectory. Based on the above discussion, the solution of the optimizing problem now centers on the following dynamic equations:

ξ(x1(t0)) )

ξ(x1(tf)) )

x2(tf) - x2(t0) )

f2(x1(t0),x2(t0),u′(t0))

(6)

f1(x1(t0),x2(t0),u′(t0)) f2(x1(tf),x2(tf),u′(tf))

(7)

f1(x1(tf),x2(tf),u′(tf)) f (x ,x ,u′)

∫xx(t(t)) ξ(x1) dx1 ) ∫xx(t(t)) f2(x1,x2,u′) dx1 1 f

1 0

1 f

1 0

1

1

2

(8)

The above three conditions can be fixed by the boundary conditions, x1(t0), x2(t0), u′(t0), x1(tf), x2(tf), and u′(tf). After determining the parameters a0′, a1′, and a2′, one can determine the desired internal state trajectory, ξ*(x1). Then we calculate u′(x1), which will drive the internal state trajectory ξ(x1) of the model to track the desired trajectory ξ*(x1) tightly. In the mean time, to map u′(x1) into u′(t), the following equation needs to be coupled with eq 4:

dt 1 ) dx1 f1(x1,x2,u′)

(9)

A conceptual flow diagram of the modified two-step method for calculating the optimizing variable, u′(t), is illuminated by Figure 1. For some practical applications of the modified two-step method described, x2(tf) may be replaced by tf. The calculation algorithm of the modified two-step method is as follows:

Ind. Eng. Chem. Res., Vol. 35, No. 7, 1996 2249

For the given problem, the following minimization problem is proposed:

objective x1f

(||u′(x1) - u′max||2 × ∑ x

min {g(x2(tf),x2(t0)) + [ x2(tf)

10

H1(u′(x1))) + (||u′min - u′(x1)||2H2(u′(x1)))]} (10) subject to eqs 4-8 with the penalty functions:

Figure 1. Computation of an optimal temperature trajectory for operating a batch polymerization reactor based on the modified two-step method.

1. Divide the interval of x1 along [x1(t0), x1(tf)] into 100 or 200 subsections, ∆x1. A constant value of the optimizing variable, u′, over each subsection, ∆x1, is assumed. 2. Integrate the system equations (eqs 4 and 5) from x1 to x1 + ∆x1 with the assumed constant value of u′. 3. Calculate the optimizing variable u′ for the next interval, ∆x1, based on the states x1 + ∆x1 and x2(x1 + ∆x1) and the desired internal trajectory, ξ*(x1+∆x1). A one-dimensional search for u′ is involved based on eq 4. 4. Go back to step 2. The whole calculation proceeds until x1 ) x1(tf). The set-point trajectory, u′(x1) and u′(t), over the course of [x1(t0), x1(tf)] can be obtained. For the problem considered, the initial conditions, x1(t0), x2(t0), and u′(t0), are usually available. However, the way to decide the boundary conditions, x1(tf), x2(tf), and u′(tf), involves different objectives to be optimized. Two cases are considered here and the possible solutions are stated as follows. Case 1. The desired performance index depends on given variables x1(tf), x2(tf), and u′(tf). For such a problem statement, the computational procedure given above (steps 1-4) can be followed accordingly unless the obtained u′(x1) gets stuck at u′min or u′max. This condition occurs at some points between x1(t0) and x1(tf). At these points, even an upper or lower bound of the optimizing variable u′(x1) cannot provide a value of ξ(x1) to meet the desired value of ξ*(x1). The way to deal with this situation is to set u′(x1) ) u′min or u′max. However, the desired internal state trajectory ξ*(x1) must be recalculated. This requirement involves recalculating the coefficients in the ξ*(x1) based on the current states, x1, x2(x1), and u′(x1), and the final conditions, x1(tf), x2(tf), and u′(tf). Case 2. The desired performance index depends on x1(tf), x2(tf), and u′(tf), among which the values of x1(tf) and u′(tf) are given but the value of x2(tf) is a free variable to be specified optimally.

H1(u′(x1)) )

{

H2(u′(x1)) )

{

0, if u′max g u′(x1) K, if u′max < u′(x1) 0, if u′(x1) g u′min K, if u′(x1) < u′min

(11)

(12)

By solving the formulated constrained problem, the optimizing variable, u′(x1), will be kept away from the constraints of u′. In this way, the operational capability of u′(x1) will be improved. Furthermore, if the obtained absolute value of the rate change of u′(t) is greater than a given value, an additional penalty function can be added in eq 10 to cope with this problem. The above minimization problem was solved by a nonlinear programming package (Ladson et al., 1989). Although the modified two-step method developed here is limited to a single-input/single-output (SISO) trajectory-tracking problem, it can be applied either offline or on-line because the computation time for each subsection is quite short. For on-line implementation, it will be advantageous to start the computational algorithm based on the initial measurement of the states of the reacting mixture, which is usually not in a planted condition. 3. Modified Globally Linearizing Control (MGLC) Structure The SISO MGLC structure (Chang and Huang, 1994) was applied to track the obtained trajectory, u′(x1) or u′(t). If the process model is exact, then the MGLC structure provides a perfect model inversion control law for trajectory tracking. A tuning method for the SISO MGLC Structure can be found in our previous works (Chang and Huang, 1994; Chang and Hseih, 1995). 4. Experimental System A schematic diagram of the experimental apparatus is provided in Figure 2. The reactor is a 0.0122-m3 (10L) stainless-steel cylindrical vessel with a 0.0254-m (1in.) drain centered at the bottom. The vessel is composed of an annular jacket with baffling and internal helical coils. The reacting mixture is blended with an agitator (diameter ) 0.22 m) at one end of which is a pitched flat-blade paddle and at the other end is a pitched anchor-type paddle. Adoption of an anchor-type paddle at the lower part of the agitator is to maintain perfect mixing. The rotational speeds of the agitator are 0-120 rpm. Detailed specifications of the reactor vessel are listed in Table 1. A Pyrex Graham condenser with a drip tip was attached to a port of the reactor head. The condenser was used to condense any vapor escaping from the reactor. Tap water (roughly 25 °C) was used as the coolant on the shell side of the condenser. A check valve

2250

Ind. Eng. Chem. Res., Vol. 35, No. 7, 1996

Figure 2. Schematic diagram of the experimental apparatus. Table 1. Specifications of Reactor Vessel DT (m) ) 0.24 D (m) ) 0.22 Djo (m) ) 0.26 dco (m) ) 2.13 × 10-2 Dc (m) ) 0.18 V (m3) ) 1.22 × 10-2

H (m) ) 0.27 N (rps) ) 7.2 × 103 Dji (m) ) 0.24 dci (m) ) 1.58 × 10-2 Nc ) 6

was installed at the outlet of the condenser to prevent any back-flow of water from filling in the vapor-trap bottle. Two vessels, V1 and V3, were used to provide two temperature levels of the coolant. Another vessel, V2, was installed to serve as the heating medium. The three vessels were about 0.08 m3. An overflow pipe was inserted to maintain a constant head for each vessel. During the reaction, the temperature of water in vessel V1 could be adjusted with an on/off solenoid valve. Meanwhile, the electric heater installed in vessel V2 was adjusted with a solid state relay (SSR) to maintain the desired temperature of the heating medium. Computer control and the installed instrumentations for the experimental apparatus are both shown in Figure 2. The computer (IBM PC/AT type) was connected to a PCL-718 A/D interface card and a PCL-726 D/A interface card (both from Advantech Corp.). These two interface cards used 12-bit converters; therefore, the digital signals were 12 bits. The analog signals from the measuring elements were amplified and conditioned by AT740 (4-20 mA/0-5V) modules (from Acel Corp.). The data acquisition software was developed by us. Heating of the reaction mixture was achieved by adjusting the flow rate of hot water, Fwj, through the jacket side of the reactor with a 3.72-kW centrifugal pump. Heat removal was carried out by manipulating the flow rate of coolant, Fwc, through the internal helical coil with another 3.72-kW centrifugal pump. Two pneumatic control valves, two flowmeters, and two rotameters were installed in the flow control loops. At the inlets and outlets of the jacket and the coil, four temperature sensors were located (0-120 °C resistance temperature detector (RTD); accuracy (0.1-0.3 °C). The reactor temperature was measured by a RTD of the same type (accuracy (0.1 °C). The volumetric flow rates

inside the jacket and the coil were measured by two Compak flow transmitters (MK508 and MK2350 from Signet Industrial, Inc.) which output a 4-20-mA signal proportional to the flow rate. 5. Dynamics of the Control Elements Compared with the dynamics of the reactor, jacket, and coil, the dynamics of the control elements (the control valve, the RTD’s, the control valve pressure transducer, and the heater) were very fast. Therefore, these dynamics were neglected in the model development of the system. The steady state input/output behavior of most of these elements showed that a linear relationship existed in all signals except the one between the air pressure and the flow for the two air-toopen control valves. This nonlinearity was accounted for in the controller system. We used the quadratic calibration equation to calculate the actual flow rate from the corresponding digital signal and vice versa. 6. Mathematical Model of the Reactor System 6.1. The Reaction System. The reaction system considered in the experiment is

maleic anhydride + n-hexyl alcohol f hexyl monoester of maleic acid This reaction is a strongly exothermic reaction (∆H ) -33.5 mJ/kmol). Although the upper temperature of this reaction is less critical than for other reactions, because of reversibility of the reaction and formation of diester, temperatures exceeding 100 °C should be avoided (Hugo, 1981). Furthermore, for the H/C system designed in Figure 2, a suitable temperature difference needs to be provided in carrying out heat transfer. Therefore, the chosen constraints on Tb,min and Tb,max are given in Table 2. 6.2. Heat-Transfer Coefficient Correlations. In developing the model based control law (MGLC structure), the overall heat-transfer coefficients, Uj(Fwj) and Uc(Fwc), are needed. However, exact correlations for

Ind. Eng. Chem. Res., Vol. 35, No. 7, 1996 2251

Figure 3. Achievement of a constant temperature with a PI controller in a CSTR operation. Table 2. Model Parameters and Operating Conditions of the Batch Reactor k (m3/(kmol‚s)) ) 1.37 × 1012 exp(-12628/T) -∆H/FCp (m3‚K/kmol) ) 16.92 CA0 (kmol/m3) ) 4.55 CB0 (kmol/m3) ) 5.34 Tb,min (°C) ) 55 Tb,max (°C) ) 75

Figure 4. Achievement of a constant temperature via a MGLC structure in the batch reactor.

estimating overall heat-transfer coefficients in a specific real batch esterification reactor are not easily available. To obtain these two correlations provided by the experimental apparatus shown in Figure 2, the following experiments were designed. (a) Determination of Uj(Fwj) from a Continuous Stirred Tank Reactor (CSTR) Operation. The experimental apparatus shown in Figure 2 was arranged as a CSTR operation. Because we cannot afford enough quantities of hexyl monoester of maleic acid, the product of the reaction, water was adopted as the process fluid in conducting the following experiments. A conventional proportional-integral (PI) controller, was adopted to control the temperature of the process fluid at a given set point. Because suitable insulation coated around the outside surface of the reactor surface, the heat loss from the reactor surface could be neglected. Therefore, the following unsteady state energy balance equation on the process fluid can be formulated:

dTb ) dt qwFwCpwTwi + UjAj(Tj - Tb) - qwFwCpwTb (13)

The correlation, Uj(Fwj), by adopting water as the process fluid, can be obtained through a CSTR operation. For a different process fluid, the heat-transfer coefficient under a similar agitation condition generally depends on the physical properties of that fluid. Based on the correlations given by Bondy and Lippa (1983), the ratio between hb,water and hb,ester can be estimated. In this way, estimation of the correlation Uj(Fwj) for monoester from the correlation for water can be obtained. (b) Determination of Uc(Fwc) from a Batch Reactor Operation. For the H/C system provided in the batch reactor (Figure 2), if the batch side is filled with hexyl monoester of maleic acid, then the energy balance equation is formulated:

FbCpbV

dTb ) UjAj(Tj - Tb) + UcAc(Tc - Tb) dt

If a given set point, Tb, is achieved by the MGLC structure, then the next equation holds:

FwCpwV

If the given set point of the process fluid is achieved, then the next equation holds:

Uj(Fwj) )

-qw(Twi - Tb)FwCpwV

(14)

Aj(Tj - Tb)

Because there is a correspondence between qw and Fwj, we can obtain the correlation Uj(Fwj) by adjusting the feed rate of water, qw. Figure 3 depicts an experimental result.

(15)

Uc(Fwc) )

-Uj(Fwj)Aj(Tj - Tb)

(16)

Ac(Tc - Tb)

On the basis of the above equation, the ratio of Uc(Fwc) over Uj(Fwj) is a function of the set point, Tb. Therefore, one can track a set of different set points by suitable paring of Fwc and Fwj decided by the control law experimentally. Then the correlation, Uc(Fwc), can be determined by the set point, Tb, and the correlation, Uj(Fwj). An experimental result is shown in Figure 4. Observation of the experimental data, nonlinear curve fitting (Coulson and Richardson, 1977), was applied in fitting the experimental data (Figure 5). The results are

2252

Ind. Eng. Chem. Res., Vol. 35, No. 7, 1996

Figure 5. Fitting the experimental data for the overall heattransfer coefficients.

1 1 ) Uj 2.673F

1 + 0.224 5.303

(17)

1 0.822

(18)

wj

1 1 ) Uc 268.8F

0.47

+

wc

6.3. Reactor Dynamics. The mathematical model was developed on the basis of the following assumptions: mixing of the reacting mixture is a perfect one; the temperatures of the fluid in either the jacket side or the coil side are represented by their mean value measured at the inlet as well as the outlet; energy dissipation induced by shaft work is neglected; heat loss to the surroundings from the surface of the reactor vessel is ignored. Consequently, species balances for the maleic anhydride, n-hexyl alcohol (Westerterp et al., 1984), and an energy balance for the reactor give the set of ordinary differential equations:

dCA ) -k(Tb)CACB dt

(19)

CB ) CB0 - (CA0 - CA)

(20)

dTb (-∆Hr)k(Tb)CACB ) + R1(Uj)(Tj - Tb) + dt FbCpb

where

UjAj FbCpbV

R2 )

UcAc FbCpbV

time and/or suitable reaction severity can be embedded in the yield optimization simultaneously. Incidentally, the constraints of the operating condition are shown in Table 2. Applicability of the proposed modified two-step method will be examined in the following case studies. The states of the reaction system are

xT ) [x1 x2] ) [CA Tb]

(22)

Model parameters and experimental conditions are shown in Table 2. 7. Calculation of an Optimizing Operating Trajectory The methodology for optimizing an operating trajectory for a batch reactor was developed in section 2. For the system under consideration, the yield optimization problem is to achieve the performance index, 99% final conversion of MA (X ) (CA0 - CA)/CA0). The least end

(23)

For this reaction system (n ) 2), Tb is chosen as the optimizing variable. Meanwhile, the internal state trajectory is

ξ(CA) )

dt 1 ) dCA -k(Tb)CACB

(24)

The following three conditions can be used in determining the parameters a0′, a1′, and a2′ if ξ(CA) is chosen as a second order polynomial.

ξ(CA0) )

1 -k(Tb0)CA0CB0

(25)

ξ(CAf) )

1 -k(Tbf)CAfCBf

(26)

tf )

R2(Uc)(Tc - Tb) (21)

R1 )

Figure 6. Optimizing trajectories Tb*(t) obtained by the modified two-step method.

∫CC

Af

A0

ξ(CA) dCA

(27)

For the defined optimization problem, two case studies and their solution procedure were discussed in section 2. Here, the details of the reaction system under consideration are discussed as follows: Case 1. The desired performance index depends on the given values of CA(tf), Tb(tf), and tf. Achievement of the performance index can be termed only as an optimizing path rather than an optimal path, as only the yield target is considered. Because an exothermic reaction is to be carried out in a batch reactor, a smaller value of tf signifies a larger intensity of reaction to be proceeded. This effect can be observed in the following two studies. Given Tb*(0) ) Tb,max, an optimizing trajectory can be obtained (line 1 in Figure 6; path 1) if tf,desired is set to be 166.7 min. As the internal state trajectory calculated on the basis of the value of tf,desired (166.7 min) is not achievable, the calculated Tb*(t) always gets stuck at the upper tem-

Ind. Eng. Chem. Res., Vol. 35, No. 7, 1996 2253 Table 3. Boundary Conditions for Determining an Optimizing Path (Path 2)

tf,desired (s) Tb* (K) CA (kmol/m3) X (%)

initial

end of stage 1

end of stage 2

end of stage 3

348 4.55 0

6000 348 0.455 90.0

6000 348 0.137 97.0

6000 348 0.045 99.0

perature constraint. Although the internal state trajectory is recalculated along the whole path, a longer time course results (tf,actual ) 239.6 min). A theoretical conclusion is that an upper bounded isothermal temperature is a minimum end-time policy for a batch reactor. The result shown in Figure 6 concludes that the proposed modified two-step method can obtain this theoretical trajectory with ease. However, if tf,desired is set to be 333.3 min, we can calculate another trajectory (line 2 in Figure 6; tf,actual ) 342.6 min). For this case, Tb* is stuck at 75 °C from t ) 0 to t ) 83.5 min. Observing the time profile of X over the interval from 0 to 50 min (lines 4 and 5 in Figure 6), the reaction is quite intensive. Let us measure the reaction severity, θ, as the reciprocal of the time required in achieving 80% conversion. For these two cases, θ is about (1/40) min-1. The heat-transfer capability determines whether we should adopt these two trajectories as the set point for our experimental apparatus. If the heat-transfer capability is not enough, it is usually recommended to adopt a larger value of tf to provide an achievable trajectory. However, due to uncertainties of the identified heat-transfer dynamics and coefficients, one needs to make experiments to obtain a suitable value of tf in reality. To reduce the time duration in which the calculated temperature trajectory will get stuck at the upper bound of temperature while the end time tf shall be kept less than that required for line 2 of Figure 6, the whole course of the reaction was divided into three stages. Moreover, the boundary conditions for the intermediate stages were adjusted by trial and error based on the time profile of X given by line 5 of Figure 6. The result is shown in Table 3. In essence, these conditions can be calculated by suitable programming. The corresponding time profiles of Tb* and X based on the conditions given in Table 3 are found in Figure 6 (lines 3 and 6). Although reaction severity is still quite intensive (θ is about (1/56) min-1), Tb* (line 3 in Figure 6) still gets stuck at the upper bounded temperature, 75 °C, for much smaller duration (0-12.8 min). Therefore, this temperature trajectory (path 2) will be adopted in the later experimental studies. Case 2. The performance index depends on given values of CA(tf) and Tb(tf) but the value of tf is a free variable to be specified optimally. Operational safety is an additional objective. Reaction severity was not embedded in optimizing an operating trajectory in the previous studies (case 1). Now this is an additional objective to be considered. In essence, the temperature control is considered as a setpoint trajectory rate-tracking problem (Bhat et al., 1991). A temperature set-point trajectory similar to the response of a first order process represents qualitatively the form of allowable rate of temperature rise suggested by Liptak (1986). The time constant of the first order response can be a tuning parameter (Bhat et al., 1991).

Therefore, we adopted the following first order trajectory as an experimental path (path 3; line 2 in Figure 7) after considering the capability of the experimental apparatus:

Tb*(t) ) 348 - 18 exp(-t/2000)

(28)

In view of the time profile of X corresponding to this first order trajectory (line 5 on Figure 7), reaction severity (θ is about (1/78) min-1) is much less than those shown in Figure 6. To show the applicability of the proposed modified two-step method after considering the objectives of reaction severity, minimum end-time operation, and the rate change of Tb*(t), the following two paths were calculated. Here, the whole course of the reaction was divided into three stages based on the same reason cited for case 1. Over each stage, Table 4 gives suitable boundary conditions by trial and error. The optimization problem (eq 10) was executed with penalty functions under constraints of Tb*min and Tb*max (eqs 11 and 12) and the absolute rate change of Tb*(t) over the first stage only. Then the time profile of Tb* (line 3 in Figure 7) can be obtained. The corresponding time profile of X is also displayed (line 5). It can be observed that Tb* gets stuck at 75 °C over stages 2 and 3. This trajectory resembles the first order trajectory described by eq 28. To be immune from this problem, one can obtain another temperature trajectory (line 1 in Figure 7; path 4) by adding all the penalty functions over stages 2 and 3. However, a much longer operation time is required than what is needed in line 3 of the same figure. This trajectory will be adopted as a path in checking the tracking capability of the MGLC structure. 8. Synthesis of the Control Law for the Experimental System The trajectory-tracking problem can be coped with by using the MGLC structure (Chang and Huang, 1994; Chang and Hsieh, 1995) if the process is a minimum phase. In this experimental study, a SISO MGLC structure was designed to prove the theoretical performance studied (Chang and Huang, 1994). If a tight trajectory tracking of reference temperature is desired in a batch reactor, both heating and cooling of the process are necessary. According to the heat balance equation (eq 13), when two manipulated variables, Uj and Uc, are used to control a single output (Tb), the control system will be excessively determined. A way of solving this difficulty is to introduce a single parametric variable u (Jutan and Uppal, 1984) such that

Uj ) (Uj,max - Uj,min)u + Uj,min Uc ) (Uc,min - Uc,max)u + Uc,max

(29)

Clearly, u ) 0 represents the maximum cooling and u ) 1 represents the maximum heating of the system. This parametric variable is used as a single control variable for this tracking problem. After substitution of eq 29 into eq 21 and rearrangement, one can obtain the following equation.

dTb ) γk(Tb)CACB + (a1 + a2Tb) + (b1 + b2Tb)u dt (30) in which

2254

Ind. Eng. Chem. Res., Vol. 35, No. 7, 1996

γ ) (-∆Hr)/FbCpb a1 ) (Uj,minAjTj + Uc,maxAcTc)/FbCpbV a2 ) -(Uj,min + Uc,maxAc)/FbCpbV b1 ) [(Uj,max - Uj,min)AjTj + (Uc,min - Uc,max)AcTc]/FbCpbV b2 ) -[(Uj,max - Uj,min)Aj + (Uc,min - Uc,max)Ac]/FbCpbV (31) Representing the process dynamics, eqs 20 and 30 in a state variable form with x ) [CA Tb] and y ) Tb, we have

f(x) )

[

-k(Tb)CACB γk(Tb)CACB + (a1 + a2Tb)

[

0 g(x) ) b + b T 1 2 b

]

]

h(x) ) Tb

(32)

Table 4. Boundary Conditions for Determining an Optimizing Path (Path 4)

(33) (34)

Because

〈dh,g〉 ) b1 + b2Tb * 0

Figure 7. Optimizing trajectories Tb*(t) obtained by the modified two-step method.

(35)

the nonlinear system under consideration is input/ output linearizable with the relative order ) 1. Then, the state feedback transformation is

tf,desireda (s) Tb* (K) CA (kmol/m3) X (%)

v ) Ω(x,u) )

∑ βkL f(h) + (-1) k)0 k

end of stage 1

end of stage 2

end of stage 3

330 4.55 0

2876 348 1.365 70.0

8179 348 0.228 95.0

9276 348 0.045 99.0

a t f,desired was determined by eqs 10-12 with the constraint on the rate change of Tb* (if dTb*/dt g 0, then (dTb*/dt)max ) 0.01 K/s; if dTb*/dt e 0, then |(dTb*/dt)|max ) 0.05 K/s).

u ) Ψ(CA,CB,Tb,v) v - β0Tb - β1[γk(Tb)CACB + (a1 + a2Tb)] )

r)1

initial

β1(b1 + b2Tb)

(36)

The resulting v-y system is described by

or equivalently

β0Tb + β1

r)1

v-



βkLkf(h)

k)0

u ) Ψ(x,v) ) (-1)

r-1

βr〈dh,adr-1f(g)〉

(37)

r)1

βk

dky dtk

)v

(38)

The parameters β0 and β1 are completely arbitrary, meaning that the v-y system can have arbitrarily placed poles and, therefore, specifies bound-input/boundoutput (BIBO) stability characteristics (Kravaris and Chung, 1987). By straight computation, the following equation can be obtained.

Lf(h) ) γk(Tb)CACB + (a1 + a2Tb)

(39)

The relation between the transformed variable v and the manipulated variable u is given by

v ) β0Tb + β1[γk(Tb)CACB + (a1 + a2Tb)] + β1(b1 + b2Tb)u (40) or equivalently

dTb )v dt

(42)

Based on the MGLC structure, v is selected as

v(t) ) Kc(Tb*(t) - Tb(t)) +

transforms the input/output system into

∑ k)0

(41)

βr〈dh,adr-1f(g)〉u

r-1

Kc τI

∫0t (Tb*(t) - Tb(t)) dτ +

β0Tb*(t) + β1

dTb* (t) (43) dt

8.1. On-Line Minimization of Energy Consumption. Based on the H/C system shown in the experimental apparatus (Figure 2), tracking a temperature trajectory by heating and cooling simultaneously, although effective, is undesirable energywise. However, observation of the heat balance equation (eq 21) and eq 29 reveals that whether the reactor system can be operated with the least energy consumption hinges on suitable choice of the terms Uj,max and Uc,max for a given manipulated variable, u (eq 29), at each sampling time. That is, if we set fixed values of Uj,max and Uc,max over the course of operation, then no improvement of energy saving can be achieved. Thus, to be energy saving, an on-line computational algorithm for finding the desired optimal values of Uj,max and Uc,max at each sampling

Ind. Eng. Chem. Res., Vol. 35, No. 7, 1996 2255

Figure 8. Controller and process for the batch reactor (SISO MGLC structure).

instance was proposed by formulating the following minimization problem.

objective

[ (|

min ω1 Uj,Uc

dTb* - γk(Tb)CACB - R1(Uj)(Tj - Tb) dt

|)

R2(Uc)(Tc - Tb)

+ ω2(|R1(Uj)(Tj - Tb)| +

]

|R2(Uc)(Tc - Tb)|) (44) subject to

h j,max U h j,min e Uj e U U h c,min e Uc e U h c,max

(45)

The term Tb* in eq 44 is determined by the modified two-step method based on the desired performance index discussed in section 7. ω1 and ω2 are the weighting factors. Minimization of the first term in the objective function (eq 44) may only result in multiple solutions of Uj and Uc. However, inclusion of the second term in the objective function is expected to track a predetermined trajectory while maintaining the least requirement for energy consumption. By the correlations (eqs 17 and 18), the variables Uj(t) and Uc(t), which minimize the objective function, will give a set of least flow rates of heating and cooling mediums, Fwj(t) and Fwc(t). In essence, for control implementation, it is desirable to adopt a lower value of the actual manipulated variables, Fwj(t) and Fwc(t), to prevent the control action bounded. Multiplying Uj(t) and Uc(t) by scaleup factors s1 and s2, which are larger than 1, is expected to provide the suitable values of Uj,max(t) and Uc,max(t) in eq 29 for control implementation. As uncertainties of the related parameters adopted in eq 44 usually exist, without suitable scale-up factors, s1 and s2, the determined flow rates, Fwj or Fwc, based on eqs 17, 18, and 29 may be too low to be regulated by the control valve if there exists model mismatch for the energy balance equation. Under the circumstances, it is suggested that the values of s1 and s2 be increased to provide a wider operational region (a larger value of Uj,max(t) and Uc,max(t)) at the expense of energy consumption. The constrained optimization problem thus formulated was

solved by the nonlinear programming routine (Kuester and Mize, 1973). 8.2. Implementation of the Control Law. Figure 8 provides the block diagram of the SISO MGLC structure. Based on this block diagram, Fortran-based source codes were set up for the main part of the controller except the C-codes so as to call input and output (I/O) interfaces. Sampling time, ∆t, was set to be 2 s. Applicability of the MGLC structure for trajectory tracking was studied by performing the following steps in executing the source codes of each block given in Figure 8. 1. Execute the state observer with the initial conditions, CA(tk-1), and the measurement, Tb(tk-1), to obtain CA(tk). 2. Calculate the weighted sum of Tb*(tk) and its derivatives dTb*(tk)/dt numerically. 3. Execute the PI control loop. 4. Compute the external input v (eq 43). 5. Minimize the objective (eq 44) to obtain Uj,max(tk) and Uc,max(tk), and then multiply the values derived by suitable scale-up factors, s1 and s2. 6. Obtain the manipulated variable u based on eq 41. 7. Calculate the overall heat-transfer coefficients, Uj(tk) and Uc(tk), through eq 29. 8. Obtain Fwj(tk) and Fwc(tk) by the correlations eqs 17 and 18, respectively. 9. Adjust Fwj(tk) and Fwc(tk) through the control valves, FCV1 and FCV2, respectively. 10. Return to step 1. If the sampling period is 2 s, the actual duration (CPU time) required to execute computational tasks (PI controller, state observer, state feedback) is about 0.05 s. A microcomputer (IBM PC/AT type) was used. This small CPU time on a modest microcomputer reflects the computational efficiency of the nonlinear control method. Therefore, applicability of the proposed method to active control in an actual setting is acceptable. 9. Experimental Procedure Reagent-grade maleic anhydride (MA) and n-hexyl alcohol were used without any further purification. A total 10 L of reacting mixture was loaded in the reactor. Based on the operating conditions given in Table 2,

2256

Ind. Eng. Chem. Res., Vol. 35, No. 7, 1996

Figure 9. Time profile of temperature in a 100-mL adiabatic reactor.

4.455 kg of MA and 5.445 kg of n-hexyl alcohol were loaded for each run. Because this esterification reaction is highly exothermic, a great deal of heat release occurs when the concentrations of the reactants are high. As the heattransfer coefficient is limited, temperature difference is another resource to be considered. Two vessels (V1 and V3) of different temperature levels were provided to cope with this problem. To track the desired temperature trajectory of the reacting medium, we used a lower temperature resource (V3) during the initial period of the reaction. In the later period of the reaction, another higher temperature resource (V1) was applied. By our estimate, the switch time at X ) 0.3 may be a suitable choice because after that moment the intensity of reaction was decreased to some extent. A solenoid valve was used to adjust the temperature of the cooling medium in vessel V1. This was accomplished by adding cold water (about 25 °C) into this vessel if the temperature of vessel V1 was higher than a specified one. In case two factorssthe ratio of the heat-transfer coefficients (hj/hc) and the temperature difference between the batch side and jacket sidesare taken into account, a suitable temperature difference between the batch side and the coil side (about 10-15 °C) is preferable. Mixing these two reactants soon caused a temperature drop. The behavior (Figure 9) was observed in a 150-mL bottle reactor which was a jacketed and insulated glass cylinder. This bottle reactor was operated adiabatically. The temperature of the reacting medium rose again because of the heat from the reaction. For the control experiment done on the experimental apparatus (Figure 2), loading 10-L reactants took about 30 s. During the reaction, temperature trajectory tracking was executed by the program automatically except the manual switching from one vessel (V3) to another (V1) manually. At the end of the reaction, we collected a sample (100 mL) and put it into the bottle reactor. The bottle reactor was initially heated to a desired temperature with a circulating heating medium provided by a constant-temperature water bath. Then the water in the jacket of the bottle reactor was pumped out just before loading the 100-mL sample from the batch reactor (Figure 2). The time profile of the temperature for the bottle reactor containing the reacting sample was recorded. For comparison, a reference sample (fully reacted material) was also tested and

Figure 10. Tracking of an optimizing trajectory Tb*(t) (path 1). Table 5. Related Parameters Used in the MGLC Structure Kc, τ I

β0, β1

sampling time

s1

s2

1, 100

1, 36

2

1.5

1.7

Table 6. Initial Temperature of V3 and V1 as Well as the Switching Time (from V3 to V1) case study

Figure

T0 (V3) (°C)

T0 (V1) (°C)

switch time (s)

path 1 path 2 path 3 path 4

12 14 15 16

53.5 53.6 43.5 44.3

62 62 53.3 48.8

900 890 810 800

recorded on the same plot. Such a method was adopted to partially assure whether the final reacting mixture achieved a 99% conversion for each run. 10. Experimental Results and Discussion Performance of the SISO MGLC structure for trajectory tracking of the paths cited above was investigated experimentally. The related parameters used in the MGLC structure are listed in Table 5. The chosen parameters, s1 ) 1.5 and s1 ) 1.7, signify that design of the control law strikes a compromise between the problem of energy consumption and uncertainties of the adopted parameter (eq 44). Table 6 gives the initial temperatures V1 and V3 as well as the switch time from V3 to V1. (a) Isothermal Trajectory (Path 1). Figure 10 depicts the set-point and reactor temperature profiles for trajectory tracking of this path. In the same plot, the average temperatures of inlet and outlet of jacket and coil, and the corresponding heating and cooling flow rates are shown. As this figure shows, the reactor was initially at 74.6 °C, which is below the set point, 75 °C. Theoretically, this isothermal path provides a global least end-time operation if the problem of reaction

Ind. Eng. Chem. Res., Vol. 35, No. 7, 1996 2257

Figure 11. Time profiles of temperature of the sample and the reference sample.

severity can be handled with the designed H/C system (Figure 2). Two different temperature levels of cooling water (V1 and V3) operated properly. Besides the superiority of the proposed control law, two other factors contributed to the successful trajectory tracking. One was the on-line manipulation of minimization of energy consumption. The other was the on-off control by the solenoid valve to roughly maintain a temperature difference between the reacting medium and the coolant in the cooling coil. Because a set of low flow rates of heating and cooling water was estimated by the minimization programming (eq 44), the desired temperature in the vessel (V1) can be achieved by adding a least quantity of supply water (25 °C). Excellent performance of the proposed MGLC structure (Chang and Huang, 1994) in handling the trajectory-tracking problem was first demonstrated experimentally here. Meanwhile, a sample was used to verify 99% conversion. Time profiles of temperature of the sample and the reference sample are shown in Figure 11. Because the difference was minor, achievement of 99% was assured. (b) Trajectory Focused on Optimizing the Yield of Reaction (Path 2). An upper bounded isothermal path is usually not desirable because the set points always get stuck at the constraint of the process capability (H/C system). Then the trajectory (line 3 in Figure 6) can be an alternative. Figure 12 depicts the experimental result. The tracking ability provided by the MGLC structure is reliable. (c) First Order Trajectory (Path 3). Here, the MGLC structure was applied to track the first order trajectory (eq 28). Because a 90 ( 2 °C temperature was maintained by the electric heater installed in the heating vessel (V2), load disturbances existed over the course of operation. However, the performance of trajectory tracking (Figure 13) was acceptable even when a temperature difference of 4.5 °C existed between Tb(0) and Tb*(0). (d) Trajectory Focused on Optimizing the Yield of Reaction, Minimal End Time, Reaction Severity, and Rate Changes of Set Points (Path 4). The trajectory under discussion is a path depicted in line 1 of Figure 7. Here the set-point policy was applied as an on-line mode where Tb*(0) was set to be the reactor temperature at t ) 0. Then the modified two-step method was applied to calculating Tb* on-line. In this way, a mismatch of temperature between Tb*(0) calcu-

Figure 12. Tracking of a first order trajectory Tb*(t) (path 2).

Figure 13. Tracking of an optimizing trajectory Tb*(t) (path 3).

lated off-line and Tb(0) can be avoided. The time profile of this trajectory as depicted in Figure 14 signifies a tough path to be tracked. However, the experimental result is acceptable. Finally, the availability of the proposed energy-saving method (eqs 44 and 45) was addressed. The reactor was loaded with the fully-reacted material. The controller

2258

Ind. Eng. Chem. Res., Vol. 35, No. 7, 1996

Figure 14. Tracking of an optimizing trajectory Tb*(t) (path 4).

Figure 16. Tracking of a first order trajectory Tb*(t) in the batch reactor loaded with fully reacted material considering the energysaving problem.

other hand, for considering the energy-saving problem (NP block was added) with s1 ) 1.1 and s2 ) 1.1, the result is shown in Figure 16. Comparison of these two figures shows that the latter used a smaller flow rate of cooling water. The difference is noticeable over the course of operation. By calculating the heat load of the coolant, heat removal was 207 kJ for Figure 16 compared with 221 kJ for Figure 15. A 6.8% energy saving can be observed. Although the quantity of the energy saving is limited, which may be attributed to the small gradient for Uc(Fwc) shown in Figure 5, the result does follow our prediction. 11. Conclusions

Figure 15. Tracking of a first order trajectory Tb*(t) in the batch reactor loaded with fully reacted material neglecting the energysaving problem.

needs to tackle the modeling-error problem (no reaction heat was released). Experimental results shown in Figures 15 and 16 are acceptable about this point. Meanwhile, for the case when the NP block (Figure 8) was dropped for neglecting the energy-saving problem, the experimental data are plotted in Figure 15. On the

An optimization and control strategy has been proposed in this work to track the temperature trajectorytracking problem in a batch reactor system. The least end time and/or suitable reaction severity can be embedded in the yield optimization simultaneously with the proposed modified two-step method in developing a set-point policy. We adopt the reliable kinetic model of the reaction system considered alone instead of the integrated dynamic model of the whole system. Because of the simplicity of the proposed two-step method, application of the proposed method in an off-line mode or an on-line mode can be applied in a real plant. The synthesis of hexyl monoester of maleic acid reaction was tested to illustrate the applicability of the proposed optimization method. In addition, the SISO MGLC control structure was proposed to tackle the optimal trajectory-tracking problem. It was developed on the basis of the global input/ output linearization technique. The model-inversion characteristic was embedded in the MGLC structure. Therefore, achievement of the desired performance

Ind. Eng. Chem. Res., Vol. 35, No. 7, 1996 2259

index, the 99% final conversion of MA, was guaranteed properly if the trajectory was tracked tightly. Furthermore, application of the SISO MGLC structure in the experimental studies was successful though a great quantity of reaction heat was released during the initial period of reaction. Model-inversion characteristic predicted from the theoretical study was demonstrated by implementation of the MGLC control structure in the experimental studies. Energy saving is another key point based on this conventional heating and cooling system designed in the experimental apparatus. We developed an on-line minimum energy-saving method for maintaining the capability of the designed H/C system during the reaction. In this way, applicability of the advanced process control techniques can be adopted accordingly in a conventional batch reactor system. Acknowledgment We thank the National Science Council (NSC 830402-E-036-007) and Dr. T. S. Lin, President of Tatung Institute of Technology, Taipei, Taiwan, ROC, for all the support in accomplishing this work. Nomenclature a1, a2 ) coefficients of eq 30 a0′, a1′, a2′ ) coefficients of eq 5 A ) heat transfer area, m2 A0 ) frequency factor of esterification reaction, m3/(kmol‚s) b1, b2 ) coefficients of eq 30 C ) molar concentration of reagent, kmol/m3 Cp ) heat capacity, kJ/(kg‚K) CV ) control valve d ) diameter of cooling coil, m D ) impeller diameter, m Dc ) centerline diameter of internal coil helix, m Dj ) diameter of annular jacket, m Dj ) outer diameter of annular jacket, m DT ) inner diameter of vessel, m E ) primary element fi|i)1-n ) function of states, eq 3 f ) vector fields that characterize the state model F ) flow rate, m3/s g ) vector fields that characterize the state model h ) state/output map; heat-transfer coefficient, kJ/(m2‚K‚s) H ) height of vessel, m H1, H2 ) penalty functions for the internal variable u′, eq 10; transfer function -∆H ) heat of reaction, MJ/kmol I ) indicator, Figure 2 k ) rate coefficient, m3/(kmol‚s) K ) positive constant value, eqs 11 and 12 Kc ) proportional gain L ) length of coil, m Lkfh ) kth order Lie derivative of h with respect to f N ) agitator speed, rev/s Nc ) number of coil turns q ) volumetric flow rate, m3/s r ) relative order of controlled output with respect to manipulated variable s1, s2 ) scaling factors S ) switch, Figure 2 t ) time, s tk ) kth sampling time, s ∆t ) sampling time T ) temperature, K or °C T ) transmitter, Figure 2 T h ) mean value of temperature measured at inlet and outlet, °C

u ) actual manipulated variable u′ ) optimizing variable U ) overall heat-transfer coefficient, kJ/(m2‚K‚s) U h ) bounded value of overall heat-transfer coefficient, kJ/ (m2‚K‚s) v ) external input of linear closed system, eq 36 V ) volume of reacting mixture, m3 V ) vessel x ) state variable x ) vector of state variables X ) conversion y ) output variable Y ) relay Greek Symbols R1, R2 ) process parameters, eq 22 β0, β1 ) tunable parameters of input/output linearized system γ ) process parameter, eq 30 θ ) reaction severity, min-1 ξ(x1) ) internal state trajectory F ) density (kg/m3) τ ) dummy variable for integration τI ) integral time constant, s Ψ ) transformation, eq 37 ω ) frequency, Hz ω1, ω2 ) weight factor, eq 44 Ω ) input-output linearizing transformation, eq 36 Acronyms A/D ) analog to digital BIBO ) bound-input/bound-output D/A ) digital to analog H/C ) heating and cooling MA ) maleic anhydride MGLC ) modified globally linearizing control PI ) proportional-integral RTD ) resistance temperature detector SISO ) single-input/single-output SSR ) solid state relay Superscript * ) optimal Subscripts 0 ) initial; inlet A ) reagent A B ) reagent B b ) reacting mixture in the batch reactor c ) cooling coil f ) final i ) inner j ) heating jacket max ) maximum min ) minimum o ) outer r ) reaction w ) water

Literature Cited Bhat, J.; Madhavan, K. P.; Chidambaram, M. Multivariable Global Input/Output Linearized Internal Model Control of a Semibatch Reactor. Ind. Eng. Chem. Res. 1991, 30, 1541-1547. Bondy, F.; Lippa, S. Heat Transfer in Agitated Vessels. Chem. Eng. Sci. 1983, April, 62-71. Chang, J. S.; Lai, J. L. Computation of Optimal Temperature Policy for Molecular Weight Control in a Batch Polymerization Reactor. Ind. Eng. Chem. Res. 1992, 31, 861-868. Chang, J. S.; Huang, K. L. Performance Study of Control Strategies for Trajectory Tracking Problem of Batch Reactors. Can. J. Chem. Eng. 1994, 72, 906-919.

2260

Ind. Eng. Chem. Res., Vol. 35, No. 7, 1996

Chang, J. S.; Chung, S. T. On-Line Set-Point Optimization and Control of a Batch Polymerization Reactor. J. Chin. Inst. Chem. Eng. 1995, 26, 213-236. Chang, J. S.; Hseih, W. Y. Optimization and Control of Semibatch Reactors. Ind. Eng. Chem. Res. 1995, 34, 545-556. Coulson, J. M.; Richardson, J. F. Chemical Engineering, 2nd ed.; Pergamon Press: Elmsford, New York, 1977; Chapter 6, p 268. Ellis, M. F.; Taylor, T. W.; Jensen, K. F. On-Line Molecular Weight Distribution Estimation and Control in Batch Polymerization. AIChE J. 1994, 40, 445-462. Hugo, P. Start-up and Operation of Exothermic Batch Process. Ger. Chem. Eng. 1981, 4, 161-173. Isidori, A. Nonlinear Control Systems An Introduction; Springer-Verlag: Berlin, 1989; Chapter 4, pp 165-172. Jutan, A.; Uppal, A. A Combined Feedforward-Feedback Servo Control Scheme for an Exothermic Batch Reactor. Ind. Eng. Chem. Res. 1984, 23, 597-602. Jutan, A.; Rodriguez, E. S., II. Application of Parametric Control Concepts to Decoupler Design and Heating Control Design for Batch Reactor. Can. J. Chem. Eng. 1987, 65, 858-866. Kravaris, C.; Chung, C. B. Nonlinear State Feedback Synthesis by Global Input/Output Linearization. AIChE J. 1987, 33, 592603. Kravaris, C.; Wright, R. A.; Carrier, J. F. Nonlinear Controllers for Trajectory Tracking in Batch Processes. Comput. Chem. Eng. 1989, 13, 73-82. Kuester, J. L.; Mize, J. H. Optimization Techniques with Fortran; McGraw-Hill, Inc.: New York, 1973; Chapter 7, pp 286-295. Ladson, L. S.; Waren, A. D.; Ratner, M. W. GRG2 User’s Guide; Department of General Business, University of Texas: Austin, TX, 1980. Lewin, D. R.; Lavie, R. Design and Implementing Trajectories in an Exothermic Batch Chemical Reactor. Ind. Eng. Chem. Res. 1990, 29, 89-96.

Liptak, B. G. Controlling and Optimizing Chemical Reactors. Chem. Eng. 1986, 26 (May), 69-81. Palanki, S.; Kravaris, C.; Wang, H. Y. Optimal Feedback Control of Batch Reactors with A State Inequality Constraint and Free Terminal Time. Chem. Eng. Sci. 1994, 49, 85-97. Ponnuswamy, S. R.; Shah, S. L. Computer Optimal Control of Batch Polymerization Reactors. Ind. Eng. Chem. Res. 1987, 26, 2229-2236. Sage, A. P.; White, C. C., III. Optimal Systems Control, 2nd ed.; Prentice-Hall, Inc.: Englewood Cliffs, NJ, 1977; Chapter 4. Soroush, M.; Kravaris, C. Nonlinear Control of a Batch Polymerization Reactor: an Experimental Study. AIChE J. 1992, 38, 1429-1448. Soroush, M.; Kravaris, C. Optimal Design and Operation of Batch Reactors. 1. Theoretical Framework. Ind. Eng. Chem. Res. 1993, 32, 866-881. Takamatsu, T.; Shioya, S.; Okada, Y. Molecular Weight Distribution Control in a Batch Polymerization. Ind. Eng. Chem. Res. 1988, 27, 93-99. Westerterp, K. R.; Van Swaaij, W. P. M.; Beenackers, A. A. C. M. Chemical Reactor Design and Operation; Wiley: New York, 1984; Chapter VI, pp 276-277.

Received for review August 28, 1995 Revised manuscript received April 1, 1996 Accepted April 2, 1996X IE950534M

X Abstract published in Advance ACS Abstracts, June 1, 1996.