Application of Nonlinear Predictive Control to a Semibatch

Ind. Eng. Chem. Res. 2004, 43, 7303-7311

7303

Application of Nonlinear Predictive Control to a Semibatch Polycondensation Reactor Galo A. Carrillo Le Roux* and Reinaldo A. Teixeira Departamento de Engenharia Quı´mica, Escola Polite´ cnica da Universidade de Sa˜ o Paulo, Avenida Prof. Luciano Gualberto Travessa 3, 380, CEP 05508-900 Sa˜ o Paulo - SP, Brazil

In this work, a model of a semibatch poly(ethylene terephthalate) reactor in which dimethyl therephthalate (DMT) and ethylene glycol (EG) react during the transesterification stage is developed and used for an application of nonlinear model predictive control (NMPC). The reactor is connected to a packed-bed column with a partial condenser. In the model, the hydrodynamic behavior of a column section is described by a correlation valid for hydrodynamic regimes under the loading point. The description of the heating of the reactor includes the division of the equipment into different compartments, each with a given thermal inertia and with a thermal resistance assumed between them. The column is considered dry at startup. NMPC is applied to the system model, and its performance is compared to that of an optimally tuned PID controller on three different temperature trajectories. The performance of NMPC is better than that of PID. The NMPC performance is studied for different modeling gaps, and it is not seriously degraded for the cases studied. This study does not intend to prove NMPC superiority over PID in the control of a semibatch polycondensation reactor. The application of NMPC illustrates the main difficulties involved in the control of the process, as the model presented here is able to describe the main features of its dynamics. 1. Introduction Polycondensation products are obtained by a reversible step growth polymerization. To obtain high-molecular-weight products, the conversion should be as high as possible. Given that the reaction is reversible, one way to favor the conversion is by removing volatile products (called the condensate) from the system. The necessity of drawing out the condensate while avoiding the entrainment of reactants implies the connection of a separation unit to the reactor. In this sense, semibatch polycondensation processes are essentially batch reactive distillation processes. In the 1980s, many works on the optimization of batch polycondensation reactors were produced. Many of these works were collected by Gupta and Kumar1 in a reference book. The interest in the field has decreased steadily since then, with rare exceptions such as the works by Sorensen and Skogestadt2 and Robertson et al.3 It is not necessary to say that experimental applications are extremely rare and, in general, are related not to control theory but to process development and kinetics. The main reasons for the limited concern in polycondensation are the relation between the added value and scale of the products, the control difficulties, and the variety of equipment designs that can be employed. Intensive commodity polycondensation products, such as poly(ethylene terephthalate) and nylon are currently mainly produced by continuous processes. Many polycondensation products are obtained in medium/small scale, with sometimes-high added value. However, the productivity that can be obtained by optimizing the dynamic behavior of polycondensation reactors does not * To whom correspondence should be addressed. Tel.: (5511) 30912246. Fax: (5511) 30912246. E-mail: [email protected].

yield a significant financial return if the scale of production is small. If the accurate control and optimization of a specific batch polycondensation system requires a complex model describing the kinetics, vaporliquid equilibrium, separation, and hydraulics, the effort needed to obtain such a model would be considerable. The financial return would not justify this effort if the scale of production were small, even if the added value of the product were high. In addition, these processes are particularly difficult to operate because of their nonlinearity and high level of interaction among manipulated variables. As the growth of polymer is affected by the amounts of reactants and condensate in the system and by the temperature of the medium, the control of these variables is essential to reach objectives such as product quality. While performing experiments in this area, our research team has experienced control difficulties that were shown to be clearly reflected in the quality of the product.4 Recently, Waschler et al.5 performed an analysis of a simplified model to extract the main physical features of the dynamics of a continuous stirred tank reactor (CSTR) in which there is vapor-liquid equilibrium. They showed that, in many conditions, such systems exhibit multiple steady states. One of these conditions is that of a kinetically controlled process with considerable differences in the boiling points of the pure components. The case studied here presents these characteristics. The analysis by Waschler et al.5 does not exactly apply to this work, as a semibatch reactor is studied (and not a CSTR); nevertheless, it shows the degree of difficulty that can be found in controlling simultaneous reaction-equilibrium systems. Some authors have studied the control problem of polycondensation reactors assuming standard column designs (i.e., with plates and reflux drum).2,6 In practice,

10.1021/ie0343317 CCC: $27.50 © 2004 American Chemical Society Published on Web 10/05/2004

7304

Ind. Eng. Chem. Res., Vol. 43, No. 23, 2004

packed-bed columns and partial condensers are more typical for polycondensation reactors because of their flexibility and their small condensate inventory. In addition, the authors do not assume that the column and reflux drum (if the latter exists) should start dry, which polycondensation reactors usually do because they are multiproduct. Because the compositions of condensate products between different batches are sufficiently dissimilar, it would be necessary to empty the column (and reflux drum) between different products. Even if the problem does not present a short-term economic impact, it does present interesting challenges. If a general and flexible solution were obtained that could be applied to different equipment and easily tailored for different products, its economic potential would be considerable. In this work, a model of a semibatch poly(ethylene terephthalate) reactor in which dimethyl therephthalate (DMT) and ethylene glycol (EG) react during the socalled transesterification stage7 is developed and used for control studies. The model corresponds to an experimental device available in our laboratory and is suited for testing control strategies. The model is able to describe the fact that the column starts dry, which introduces discrete switches between operating modes that describe a single phase behavior or vapor-liquid equilibrium, that is, vapor withdrawal from the system. This characteristic introduces discontinuities into the behavior of some variables such as the vapor flux. This complexity, in turn, introduces additional difficulties for the controller, which should act as for disturbance rejection. A nonlinear model predictive controller (NMPC) is applied to the system model to evaluate its potential use. It is based on a very detailed model. The NMPC would be expected to have very good performance because it is able to deal with nonlinearities and it is suited for predicting discontinuities in the process variables, which are deterministic because they arise from operating conditions and not from unmeasured disturbances. The NMPC implemented serves as a benchmark that is useful to quantify the performance of any other control scheme in comparison. The implementation of the NMPC on the model provides the opportunity of studying robustness issues, as it enables the intentional introduction of modeling errors. This is important for linking the ability of the NMPC to deal with nonlinearities to robustness. In this work, an optimally tuned PID controller is compared to the NMPC in terms of performance. This study does not intend to prove superiority of the NMPC over the PID in the control of semibatch polycondensation reactors. The application of the NMPC illustrates the main difficulties involved in the control of the process, showing that the model presented here is able to describe the main features of its dynamics. The application of the NMPC to the available experimental device is not yet intended, as it would be of little practical use because the model is specific for the poly(ethylene terephthalate) system and the experimental reactor described herein. Too much effort would be necessary to implement the NMPC, and it is more important to concentrate on the development of a more general control solution.

Figure 1. Schematic representation of the experimental device.

2. System Studied The model developed and studied in this work is suited to describe an experimental reactor available in our laboratory, which is represented in Figure 1. The reactor is heated by means of an electric resistance, soldered in aluminum, at the bottom of the reactor. The heating power is manipulated remotely through a programmable logical controller (Allen-Bradley, model SLC/5). Electrical resistance is used, instead of a liquid transfer medium, because the reactor is able to operate at temperatures up to 350 °C. The reactor is equipped with a mixer with controlled speed, remote speed set point, and torque measurement (CAT, model R 100C). A packed-bed column is attached to the vapor outlet. The top of the column is connected to a partial condenser made up of a jacket through which cooling water flows and inside which vapor condenses in a bundle of vertical tubes. A pump (ISMATEC, model Reglo-Z) controls the flow of cooling water and is manipulated remotely. A total condenser is placed after the partial condenser to condense the vapor that leaves the system. Platinum resistance thermometers (RTDs) are used to perform temperature measurement throughout the system. There is one RTD inside the reacting medium and another at the partial condenser vapor outlet. The process is supervised by the Intellution FIX HMI software. 3. System Model The kinetics of the system has been described by Shin et al.8 and corresponds to the prepolymerization step of poly(ethylene terephathalate) production. No secondary reaction is assumed to take place, as our focus is on the macroscopic effects of process variables. The reactants are dimethyl terephthalate (DMT) and ethylene glycol (EG). The condensate product is methanol (M). We consider two kinds of terminal groups, Eg (glycol) and Em (methyl) and a glycol ester bridge, which regroups an ethylene glycol molecule bonded through two ester bonds to benzene rings (Z). The reaction

Ind. Eng. Chem. Res., Vol. 43, No. 23, 2004 7305

mechanism is made up of three types of reactions

Ester Interchange k1

Em + EG {\ } Eg + M k /K 1

1

Transesterification k2

Em + Eg {\ }Z+M k /K 2

2

Polycondensation k3

2Eg {\ } Z + EG k /K 3

Figure 2. Compartments assumed to model the heating behavior of the reactor.

3

Assuming a reactor volume of V, material balances were carried out for each of the species EG and M, for the terminal groups Eg and Em, for the bonds Z, and for the volatile species considered EG and M

[ (

)

[Eg][M] dEm ) V -k1 2[Em][EG] dt K1

(

k2 [Em][Eg] -

[(

)

dEg [Eg][M] ) V k1 2[Em][EG] dt K1

(

k2 [Em][Eg] -

)

( ) (

)]

(1)

)]

(2)

)]

(3)

)]

(4)

2[Z][M] K2

2[Z][M] 4[Z][EG] - 2k3 [Eg]2 K2 K3

[(

2[Z][M] dZ + ) V k2 [Em][Eg] dt K2 k3 [Eg]2 dEG ) -yEGFV + xEG,NLN + dt

[ (

V -k1 2[Em][EG] -

(

)

[Eg][M] + K1

2k3 [Eg]2 -

4[Z][EG] K3

dTr

)

dt

dM ) -yMFV + xM,NLN + dt

[(

4[Z][EG] K3

When there is no vapor-liquid equilibrium in the reactor, these flows are both set to zero. Raoult’s law is applied to predict partial vapor pressures. The onset of vapor-liquid equilibrium is detected by examining the sum of the partial pressures of the volatile species and comparing it to the total pressure. After this condition has been achieved, an equation that states that the sum of the partial pressures is equal to the total pressure is incorporated into the system of differential equations to describe the vapor-liquid equilibrium. This equation, together with the energy balance, allows for the evaluation of the molar flow rate of vapor that leaves the reactor, FV. When vapor-liquid equilibrium occurs, the temperature cannot change freely, because the equilibrium condition must be respected at any time. For instance, if there were only one volatile component when vapor-liquid equilibrium was established, the temperature would be fixed at the boiling point of the volatile component. The onset of vapor-liquid equilibrium introduces a switch in operating behavior that could seem to be produced by an unmeasured disturbance, but that is, in fact, predictable. The description of the heating of the reactor consists of dividing the equipment into different compartments, each with a given thermal inertia and with a thermal resistance assumed between them, as described in Figure 2. The energy balance for the reacting medium takes the form

)

[Eg][M] V k1 2[Em][EG] + K1

(

k2 [Em][Eg] -

1 mmrCp,mr

nr

∆Hr ri) ∑ i)1

[(Q4 - Q5 - Q11) - (V

i

FV(HVR - hVR) + LN(HLN - hLN)] (6)

)]

2[Z][M] K2

where

(5)

The factors 2 and 4 appearing in the preceding equations are due to the number of distinct ways a component can react. The components that can react in two distinct ways are EG, because it has two alcohol end groups, and Z, because it regroups two ester bonds.8 The reaction (k1, k2, and k3) and equilibrium (K1, K2, and K3) constants were extracted from Shin et al.8 FV is the molar flow rate of vapor that leaves the reactor, and LN is the liquid flow rate from the final column section, just above the reactor, that flows back to it.

vap HVR ) yEG[∆Hvap EG + Cp,EG(Tr - Tref)] + yM[∆HM + Cp,M(Tr - Tref)] (7)

hVR ) yEGMWEGCp,mr(Tr - Tref) + yMMWMCp,mr(Tr - Tref) (8) HLN ) xEG,NMWEGCp,EG,N(TN - Tref) + xM,NMWMCp,M,N(TN - Tref) (9) hLN ) xEG,NMWEGCp,mr(Tr - Tref) + xMMWMCp,mr(Tr - Tref) (10)

7306


∆Hvap i , MWi, and Cp,i are the heat of vaporization, the molecular weight, and the specific heat of component i, where i ) EG and M. ∆Hri represents the heat of reaction ri. Q4, Q5, and Q11 represent the heat exchanged by the reacting medium with the bottom, the side, and the cover of the reactor, respectively, following Figure 2. The modeling of the fixed-bed column assumes that it is dry in the beginning and is based on the MESH equations. This is possible because the column is divided into a finite number of sections, with a given height equivalent to a theoretical plate (HETP).9 The onset of vapor-liquid equilibrium is detected in each stage. The hydrodynamic behavior of a column section is described by Billet’s correlation,10 valid for hydrodynamic regimes under the loading point. This equation asserts a relation between the liquid holdup and the liquid flow rate

Lj ) βUj3

(11)

where

where β is evaluated through the expression

β)

MWm,j2g 12ηLj a2Za3A2FLj (ah/a)2

(j ) 1, 2, ..., N) (12)

(13)

where Vj and Lj are the rates of vapor and liquid flow

dEGj ) xEG,j-1Lj-1 + yEG,j+1Vj+1 - xEG,jLj - yEG,jVj dt (j ) 1, 2, ..., N) (14) dMj ) xM,j-1Lj-1 + yM,j+1Vj+1 - xM,jLj - yM,jVj dt (j ) 1, 2, ..., N) (15) of EG and M leaving each column section j. xEG, xM and yEG, yM are the molar fractions of EG and M in the liquid and vapor phases, respectively. As in the reactor, the onset of vapor-liquid equilibrium is tested in each of the sections. The energy balance for each column section j is given by

dTj 1 (V H1 + ) dt (Mr,jCp,r,j + EGjCp,EG,j + MjCp,M,j) j+1 j+1 2 Lj-1Hj-1 - VjH3j + Qc,j) (j ) 1, 2, ..., N) (16)

and

Qc,j

{

) 0 if j > 1 * 0 if j ) 1

Qc,j ) UgA(Tj - Tcp)

(17)

1 Hj+1 ) yEG,j+1[∆Hvap EG + Cp,EG(Tj+1 - Tj)] +

with ηLj , FLj , and MWm,j representing the viscosity, density, and molecular weight, respectively, of the liquid mixture in section j of the fixed-bed column. g is the acceleration of gravity, a is the specific surface area of the packing, ah is its hydraulic area, and Za the height of section j. The ratio ah/a is assumed to be 1. A is the cross-sectional area. The global and partial (for components EG and M) material balance equations, at each column section, are given by

dUj ) Lj-1 + Vj+1 - Lj - Vj (j ) 1, 2, ..., N) dt

Figure 3. Tr and Tc as functions of time for fixed Qr and qw.

yM,j+1[∆Hvap M + Cp,M(Tj+1 - Tj)] (18) 2 Hj-1 ) xEG,j-1Cp,EG(Tj-1 - Tj) + xM,j-1Cp,M(Tj-1 - Tj) (19) vap H3j ) yEG,j∆Hvap EG + yM,j∆HM

(20)

In the energy balance for the partial condenser, the heat transfer with the cooling water is considered. The partial condenser is modeled as a column section, and an additional energy balance is performed in the cooling water jacket, given by

Qc,j dTcp qw(Tw - Tcp) + ) dt Vcp FwaterVcpCp,water

(21)

where Fwater, Vcp, Cp,water, Tw, and qw are, respectively, the density, volume, specific heat, temperature, and flow of the cooling water. The manipulated variables are Qr, the electrical power introduced to the resistance at the bottom of the reactor, called u1, and qw, the cooling water flow through the partial condenser jacket, called u2. The controlled variables are the reactor temperature (Tr) and the condenser temperature (Tc). It is possible to infer the composition of the vapor leaving the system from the condenser temperature. Thus, the condenser temperature set point (in the transesterification stage) should be close to the methanol boiling point, to avoid ethylene glycol losses. Trajectories for Tr and Tc define portable operating policies. Their application depends on adequate tracking control implementation. Henceforth, the following parameters for the system are assumed: the initial amount of DMT is 2 mol, the initial amount of EG is 4 mol, the number of sections in the column is 5 (the partial condenser is numbered 1), the initial temperature (all compartments) is 413 K, and the air temperature is 298 K. The model was implemented in Matlab and solved by the ode15s routine. In Figure 3, Tr and Tc trajectories


Figure 4. Temperatures in the sections of the column as a function of time for fixed Qr and qw.

Figure 6. Tr as a function of time for Qr ) 330 and 430 W.

Figure 5. Vapor flow in the sections of the column as a function of time for fixed Qr and qw.

are presented for fixed Qr (330 W) and qw (40 mL/min). It can be seen that, because Tc equals the methanol boiling point (335.6 K), no EG leaves the system during operation for these conditions. In Figure 4, the temperature profile along the column as a function of time can be seen. In Figure 5, the vapor flows from the different sections of the column are presented as a function of time. Although there are apparent discontinuities in vapor flow (Figure 5), no unmeasured disturbances are taking place in the system. V5 is the vapor flow from the first column section. It rises very fast as soon as the liquid in the section reaches the bubble-point temperature (Figure 4). In simple terms, this means that the liquid begins to boil. When the vapor leaves section 5, it reaches section 4, which begins to be filled with the vapor that condenses there because of the low initial temperature of the section. When the liquid in section 4 reaches its bubble-point temperature (Figure 4), vapor rises from this section (V4) and begins filling section 3, and so on, until section 1 (the partial condenser) is filled by liquid. Because the bubble-point temperature depends on composition, Figure 4 illustrates the fact that the composition in the column becomes enriched in methanol as we approach the partial condenser (section 1). Vapor does not leave the system as long as the temperature in section 1 (partial condenser) does not reach the methanol boiling point. Because the rate of reaction decreases in time because of the reactant

Figure 7. Tc as a function of time for Qr ) 330 and 430 W.

consumption, methanol withdrawal decreases and tends asymptotically to zero. The time behaviors of Tr and Tc are presented in Figures 6 and 7, respectively, for two different fixed values of Qr: 330 and 430 W. It can be seen that, as the heating power increases, vapor leaves the system sooner and its composition no longer corresponds pure methanol. This can be inferred because Tc is higher than the methanol boiling point. The control of qw is critical because, if it were too low, the temperature of the partial condenser would rise, thus allowing ethylene glycol to leave the system. If the amount of EG in the system is not controlled, the quality of the product cannot be known a priori, because ethylene glycol is a reactant. The function of qw is regulation of the reflux ratio, because it sets the ratio between the vapor and liquid flows inside the column. The behavior of Tr and Tc as a function of time is presented in Figures 8 and 9, respectively, for two different fixed values of qw: 2 and 40 mL/min. It can be noticed that the two manipulated variables have similar effects on both of the controlled variables. For the control of the composition of the vapor leaving the system, the key variable is the ratio between the vapor and liquid flows along the column, which corresponds to the reflux ratio and defines the operating line of the column.

7308

Ind. Eng. Chem. Res., Vol. 43, No. 23, 2004 Table 1. Tr Trajectories Studied

Figure 8. Tr as a function of time for qw ) 2 and 40 mL/min.

trajectory 1

trajectory 2

trajectory 3

from 413 to 453 K in 10 min from 453 to 473 K in 80 min hold at 473 K until t ) 120 min

from 413 to 453 K in 20 min from 453 to 463 K in 70 min hold at 463 K until t ) 120 min

from 413 to 453 K in 30 min hold at 453 K until t ) 120 min -

output variables, y(k + j|k) represents the predicted output for time (k + j) of the process at instant k. The control horizon M is 2; the prediction horizon P is 4; and dˆ (k) represents the estimated disturbances, calculated as the difference between the output of the process and the predicted output obtained from the nonlinear reference model. The matrix Q equals equality matrix (the two variables have the same weight), and S is assumed null. The optimization problems generated by the NMPC are solved at each control interval by the fmincon routine in Matlab. This routine implements the SQP (sequential quadratic programming) method. Constraints for the manipulated variables are defined as follows

0 e u1 e 1000 0 e u2 e 100

Figure 9. Tc as a function of time for qw) 2 and 40 mL/min.

4. Control NMPC and PID controllers performance was compared on three different temperature trajectories. The trajectories are described in Table 1. As you move from trajectory 1 to trajectory 3, the slope of the temperature ramp and the final temperature become less severe. The set point of Tc (340 K) is a little higher than the pure methanol boiling point (335.6 K). This is a typical operating practice used to avoid problems that would be caused by a small temperature sensor bias. The NMPC controller implemented follows that of Henson.11 The NMPC objective function formulation has the form

min

u(k|k),u(k+1|k),...,u(k+M-1|k)

J)

To avoid numerical problems, the control of the condenser temperature is not made effective as long as the condenser temperature does not approach its set point. This is because vapor does not leave the system as long as proper conditions in the reactor and the column have not been met. As u2 (qw) would have no effect on the two controlled variables thus far, the optimization problem would be singular. A digital PID algorithm was employed to compare the control performance. The PID tuning was obtained by optimizing its parameters to minimize the sum of the integral of square errors (ISE) for both controlled variables. The initial value of manipulated variable u1 was set to zero. The control interval is 0.1 min, and the pairing of manipulated and controlled variables is Tru1 and Tc-u2. The solution for the integral constant of loop 1 would be zero if there were no constraint on the integral value. This solution leads to the minimum ISE value, but it makes the system very sensitive to disturbances because there is no integral action, and control performance would present an offset. Thus, a minimum value of 23 was established for the integral constant, and the solutions obtained for the tuning parameters of the loops are

loop 1: Kp ) 100, Ki )23, Kd ) 100 loop 2: Kp ) 69.7, Ki ) 57.4, Kd ) 9.1

T

[y(k + P|k) - ysp(k) + dˆ (k)] Q[y(k + P|k) ysp(k) + dˆ (k)] + P-1

[y(k + j|k) - ysp(k) + dˆ (k)]TQ[y(k + j|k) ∑ j)0 ysp(k) + dˆ (k)] + ∆uT(k + j|k)S∆u(k + j|k)

(22)

where ysp(k) represents the set-point trajectories for the

4.1. Base Case. In Figures 10-12, the performance of the PID and NMPC for Tr trajectories 1-3, respectively, is presented. These figures are magnifications of the original trajectories with time spans of only 40 min whereas the range of the original trajectories is 120 min. This is to emphasize four main events, where the system behavior presents difficulties in tracking the Tr trajectory: (i) The first event is at the beginning of operation, when the system must leave its initial state and follow the initial ramp. The NMPC response is quite fast, the


Figure 10. Tr trajectory 1 with PID and NMPC controllers.

Figure 13. Tc trajectory 1 with PID and NMPC controllers.


Figure 14. Tc trajectory 2 with PID and NMPC controllers.


whereas PID controller takes 5-10 min to track the trajectory (Tr-I). (ii) At about 8 min for trajectory 1, 13 min for trajectory 2, and 15 min for trajectory 3, the reactor reaches vapor-liquid equilibrium, and the behavior is similar to what would be introduced by an unmeasured step disturbance (Tr-II). When the NMPC predicts the discontinuity in operating behavior, it

promptly recovers from the situation, whereas it takes longer for the PID controller. (iii) The third event (TrIII) corresponds to the interaction of Tc stabilization and will be better understood in the analysis of the Tc trajectories. In Figure 10, a small overshoot can be identified at about 12 min, in Figure 12 at 23 min, and in Figure 11 at 20 min (in fact, in this figure the effect is not remarkable). (iv) Finally, the fourth event (TrIV) corresponds to the change in the slope of the temperature ramp. Controllers manage the discontinuity on the derivative of the trajectory tracked, at time 10 min for trajectory 1, 20 min for trajectory 2, and 30 min for trajectory 3. In Figures 13-15, the behavior of Tc when PID and NMPC controllers are employed to track trajectories 1-3, respectively, is presented. The behavior of Tc is characterized by four different behaviors: (i) a region of constant temperature, in which there is no vapor leaving the reactor (Tc-I); (ii) a region with a steep slope (Tc-II) that corresponds to the thermal transient of the column, when the reactor enters vapor-liquid equilibrium and the column begins to grow warmer until condensate begins to leave the partial condenser; (iii) a region of constant temperature of 333 K (Tc-III), during which vapor first leaves the system as it reaches the partial condenser at the methanol boiling point; and (iv)

7310


Figure 15. Tc trajectory 3 with PID and NMPC controllers. Table 2. Integral of the Square Error for Each of the Controlled Variables for the NMPC and the PID Controller NMPC trajectory 1 trajectory 2 trajectory 3

Figure 16. Tr trajectory for the cases in which the kinetic constant of reaction 1 of the reference model is changed by +20% and -20%.

PID

ISE1

ISE2

ISE1

ISE2

4.4 1.0 3.9

62.7 97.2 133.2

65.0 31.2 23.3

58.6 103.6 149.8

a region (Tc-IV) in which the condensate temperature follows its set point (340 K). The flow of cooling water interacts with the control of Tr. This corresponds to the event named Tr-III. As vapor reaches the partial condenser, it first heats the water contained in the jacket. This corresponds to the stage where condensate leaves the system at 335.6 K. Water flow to the jacket is null, as any positive flow would have the effect of lowering the condensate temperature. When water begins to flow (just before Tc-IV, when the condensate leaves the system at 340 K), its effect is to increase the flow of liquid back to the reactor. The NMPC overreacts, which is why there are small overshoots at times 12 and 23 min in Figures 10 and 12, respectively. An important effect is that the time history interferes with the amount of volatile species present in the reactor. This is the reason the onset of vapor-liquid equilibrium takes places at different times for PID and NMPC controllers, as can be observed in Figures 1315. In Table 2, the performances are compared quantitatively in terms of the ISEs for manipulated variables 1 (Tr) and 2 (Tc). The NMPC performance is better than that of the PID. 4.2. Robustness. In this section, the NMPC performance is studied for different modeling gaps. The system is unchanged, but the reference model is changed. The following parameters were varied with respect to the system model: (i) the number of sections in the column for the reference model was changed to 2; (ii) some initial temperature were perturbed, with the initial temperature of the resistance compartment being raised by 50 K, the temperature of the aluminum by 30 K, and those of the base and the side by 10 K; (iii) the mass of the copper resistance was changed by +30%; (iv) the heat transfer rate in the partial condenser, UA, was changed by +30%; and (v) the kinetic constant of reaction 1 was changed by +20% and -20%.

Figure 17. Tc trajectory for the cases in which the kinetic constant of reaction 1 of the reference model is changed by +20% and -20%. Table 3. Integral of the Square Error for Each of the Controlled Variables for Different NMPC Robustness Study base case number of sections ) 2 perturbation initial temperatures mass of copper resistance UA partial condenser kinetic constant (k1) +20% kinetic constant (k1) -20%

ISE1

ISE2

4.4 7.9 8.8 12.3 6.0 8.7 4.7

62.7 285.4 52.6 64.2 61.7 69.8 62.4

The results are presented in Table 3, and Tr and Tc trajectories are presented in Figures 16 and 17 for the case in which the kinetic constant of reaction 1 is changed by +20% and -20%. The NMPC performance is not seriously degraded for the cases studied. A simple feedback introduced in the NMPC formulation was able to reject the modeling gaps.11 The worst performance was for the gap in the number of sections of the column. 4.3. Robustness and Trajectories. In Table 4, results for the three trajectories, for the cases where the number of sections in the column is changed to 2 and the initial temperatures are perturbed, are pre-

Ind. Eng. Chem. Res., Vol. 43, No. 23, 2004 7311 Table 4. ISE for the Three Trajectories for the Cases in Which the Number of Sections in the Column Is Changed to 2 and the Initial Temperatures Are Perturbed ISE1 trajectory 1 trajectory 2 trajectory 3

number of sections ) 2 7.9 11.9 6.4

initial temperatures perturbed trajectory 1 8.8 trajectory 2 29.8 trajectory 3 197.8

ISE2 285.4 370.9 458.5 52.6 153.0 191.3

sented. It can be concluded that, as trajectories go from 1 to 3, the error in Tc increases. This is due to the fact that the production of methanol decreases and, as a consequence, so does the rate of vapor from the partial condenser. This makes the time required for the vapor leaving the system to go from the methanol boiling point to 340 K longer, so that the contribution of region TcIII (as in Figures 13-15)to the total ISE becomes more important. In addition, as the rate of vapor withdrawal from the partial condenser decreases, qw also decreases, and the optimization problem becomes more sensitive to small variations of this variable. 5. Conclusions The results obtained showed that NMPC performance is better than that of the PID. The NMPC implemented was tested for robustness, and its performance was not seriously degraded. However, this study does not intend to prove NMPC superiority over PID in the control of a semibatch polycondensation reactor. The application of the NMPC illustrates the main difficulties involved in the control of the process, as the model presented here is able to describe the main features of its dynamics. The NMPC implemented here is specific for the prepolymerization step of PET production and the experimental reactor studied. Research should be focused on the development of a simpler reference model that could be readily adapted for different equipment, products, and scales of production. An important field of research is the study of the feasibility of trajectories. Few works have been published on the conditions for a trajectory to be feasible,

on a specific equipment design, or on any theoretical condition. This would lead to the proposition of trajectories known a priori to be feasible, which would be easier to track, even by simple controllers. Acknowledgment The authors acknowledge the Conselho Nacional de Desenvolvimento Cientı´fico e Tecnolo´gico, CNPq, and the Fundaç aõ de A ˆ mparo a` Pesquisa do Estado de Saõ Paulo, FAPESP, for their financial support. Literature Cited (1) Gupta, S. K.; Kumar, A. Reaction Engineering of Step Growth Polymerization; Plenum Press, New York, 1987. (2) Sorensen, E.; Skogestad, S. Control Strategies for Reactive Batch Distillation. J. Process Control 1994, 4, 205-217. (3) Robertson, D. G.; Russell, S. A.; Lee, J. H.; Ogunnaike, B. A. Modeling and control of a batch condensation reactor. In Proceedings of the American Control Conference; IEEE Press: Piscataway, NJ, 1995; pp 1746-1750. (4) Lo´pes Noriega, G. C. Otimizaç aõ Experimental da Sıńtese de um Oligoˆmero E Ä ster. M.Sc. Thesis, Escola Politećnica da USP, Saõ Paulo, Brazil, 2001. (5) Waschler, R.; Pushpavanam, S.; Kienle, A. Multiple steady states in two-phase reactors under boiling conditions. Chem. Eng. Sci. 2003, 58 (6), 2203-2214. (6) Balasubramhaya, L. S.; Doyle, F. J., III. Nonlinear modelbased control of a batch reactive distillation column. J. Process Control 2000, 10, 209-201. (7) Ravindranath, K.; Mashelkar, R. A Polyethylene TerephthalatesI. Chemistry, thermodynamics and transport properties. Chem. Eng. Sci. 1986, 41 (9), 2197-2214. (8) Shin, J.; Lee, Y.; Park, S. Optimization of the pre-polymerization step of the polyethylene terephthalate (PET) production in a semibatch reactor. Chem. Eng. J. 1999, 75, 47-55. (9) Salimi, F.; Depeyre, D. Comparison between dynamic behaviour of a batch packed and plate column. Comput. Chem. Eng. 1998, 22 (3), 343-349. (10) Billet, R. Packed Towers in Processing and Environmental Technology; VCH Publishers: New York, 1995. (11) Henson, M. A. Nonlinear model predictive control: Current status and future directions. Comput. Chem. Eng. 1998, 23, 187202.

Received for review December 19, 2003 Revised manuscript received July 30, 2004 Accepted August 3, 2004 IE0343317

Application of Nonlinear Predictive Control to a Semibatch

Recommend Documents