Nonlinear Multirate Model-Algorithmic Control. 2. Experimental

4064

Ind. Eng. Chem. Res. 2002, 41, 4064-4074

Nonlinear Multirate Model-Algorithmic Control. 2. Experimental Application to a Polymerization Reactor Michael P. Niemiec,† Costas Kravaris,*,‡ and Robert Zand§ Department of Chemical Engineering, Macromolecular Science and Engineering Center, The University of Michigan, Ann Arbor, Michigan 48109

This paper experimentally implements the nonlinear multirate model-algorithmic controller to regulate the reactor temperature, which is measured at a fast sampling rate, and the numberaverage molecular weight of the polymer, which is measured at an extremely slow sampling rate. A detailed model of the continuous free-radical methyl methacrylate polymerization reactor was developed from first principles to construct the multirate controller algorithm. The performance of the controller was examined by varying the tunable controller parameters, where it was shown that the controller can adequately track step changes in the output set points under active input constraints. Decreasing the tunable controller parameters deteriorated the control of the number-average molecular weight but improved the control of the temperature. Overall, the nonlinear multirate controller was an efficient and effective method to deal with the difficult problem of multirate output measurements. 1. Introduction In many chemical processes, output measurements are often unavailable at the same sampling rate. This is especially common in polymerization reactors, where many important process variables related to product quality can only be measured at extremely slow sampling rates with considerable time delays. For example, the molecular weight distribution (MWD) is an important polymer property that directly influences several end-use characteristics of the polymer such as the impact strength, hardness, tack, and elastic modulus.1 In addition, the melting point and flow properties of the polymer, which determine the method and equipment for final processing, are also affected by the molecular weight of the polymer.2 Therefore, there is strong motivation to achieve accurate control of the molecular weight during the operation of polymerization processes. Because the MWD is most often obtained with a gel permeation chromatograph (GPC) that has a cycle time of 15-30 min, this is a difficult and challenging control problem. Over the last few decades, most research efforts on polymerization reactor monitoring and control have been centered on the use of the multirate extended Kalman filter (EKF), which incorporates the slowly sampled measurements into the estimation algorithm.2-5 The estimated variables are then used in a variety of control schemes to regulate the outputs. A proportionalintegral-derivative (PID) controller was used in ref 2, a nonlinear reference controller in ref 3, and model predictive controllers in refs 4 and 5. All of the aforementioned works controlled the polymer molecular weight and, in particular,3-5 involved experimental application to different types of polymerization reactors. * To whom all correspondence should be addressed. Tel: +30-610996339. E-mail: [email protected]. † Present address: Honeywell, 16404 N. Black Canyon Highway, Phoenix, AZ 85053. E-mail: michael.niemiec@ honeywell.com. ‡ Present address: Department of Chemical Engineering, University of Patras, GR-26500 Patras, Greece. § E-mail: [email protected].

Recently, several new approaches have been developed that do not involve the use of a multirate EKF. A multirate version of linear model-algorithmic control (MAC) was developed in ref 6 and applied to regulate species concentrations and pressure in a copolymerization reactor. With this method, the slowly sampled measurements are incorporated into the output predictions when they become available. In ref 7, a multirate nonlinear state estimator was presented which can directly use a nonlinear process model in the observer design without any linear approximation. Finally, in part 1 of this series of papers,8 a general multirate nonlinear model-algorithmic controller was developed. It incorporates the latest available measurements into its nonlinear structure, which both gives excellent control of the slowly sampled outputs and effectively rejects disturbances in the outputs with faster sampling rates. In part 2 of this series of papers, the performance of the nonlinear multirate MAC algorithm8 is tested through its experimental implementation on a continuous stirred tank polymerization reactor. The reaction involves the solution polymerization of methyl methacrylate (MMA) in toluene, which is initiated by azobisisobutyronitrile (AIBN). The reactor temperature and number-average molecular weight (NAMW) of the polymer are controlled by manipulating the net rate of heating or cooling to the reactor jacket and the flow rate of the inlet initiator stream. The NAMW measurement is obtained by a GPC and has a significantly longer sampling period and delay than the reactor temperature measurement. The paper starts with a description of the experimental reactor system designed and built for the implementation of the nonlinear multirate model-algorithmic controller. A mathematical model of the system is then developed. Finally, the performance of the controller is investigated in terms of steady-state regulation and tracking of step changes in the output set points. 2. Experimental System As depicted in Figure 1, the reactor is a 3 L jacketed glass vessel, which is mixed by a multipaddle agitator

10.1021/ie010560m CCC: $22.00 © 2002 American Chemical Society Published on Web 07/16/2002

Ind. Eng. Chem. Res., Vol. 41, No. 16, 2002 4065

Figure 2. Free-radical reaction mechanism.

Figure 1. Experimental system.

connected to a constant rpm motor. The heating/cooling system of the reactor consists of electrical heaters, a circulation pump, circulation piping, an electric control valve, and temperature sensors (two 0-100 °C resistance temperature detectors (RTDs)). The reactor temperature is measured by a RTD of the same type. These RTDs are connected to a National Instruments (NI) SC2042-RTD external board and further to a NI PCI-1200 board inside the computer. With this scheme, the inlet flow rate of cooling water (Fcw) and the power supplied to the electrical heaters (P) are manipulated inputs. The inlet coolant flow rate is adjusted by the electric control valve in proportion to a 4-20 mA signal sent from the computer (a NI AT-AO-6 analog output board is used for all analog 4-20 mA signals). The heater power is adjusted by a solid-state relay (SSR) which is connected to a pulse control module driver (both from Omega Engineering, Inc.). The driver module allows the simple conversion of the on/off SSR to a proportional power regulator. Thus, the average power to the heater is proportional to the 4-20 mA signal sent to the driver module. The electrical heaters that receive the power are four 1000 W Fire Rod resistance heaters in series. In the experimental system, there are two inlet streams and one outlet stream. The first inlet stream is the monomer stream, which is fed by a metering pump (model RHSY1CTC from Fluid Metering, Inc. (FMI)) at the desirable constant flow rate of 2.78 × 10-7 m3/s. The second inlet stream is the initiator stream. The initiator flow rate (Fi) is a manipulated input of the system and is adjusted by a 4-20 mA signal sent to the metering pump (model QVG50HOCTCLF from FMI). Both the monomer and initiator stream feed bottles, as well as the reactor solution, are purged of oxygen (a reaction inhibitor) by bubbling nitrogen through the solutions and agitating each with stirrers. The outlet (product) stream of the system is pumped out by a metering pump (model RHVOCKC from FMI) driven by a 4-20 mA signal. The variable-speed pump maintains the volume of the reacting mixture at a constant value of 1.20 × 10-3 m3.

Additionally, a small fraction of the product stream is diverted for use in the Waters GPC system. The GPC measures the average molecular weights of the polymer being produced online. The unit includes a Waters 600 controller, 410 differential refractometer, 610 fluid unit (pump), and a Rheodyne 10 port valve, which are controlled by the Millenium 32 software supplied with the system. After the product stream is diverted, it is diluted to a concentration of approximately 0.1% with additional solvent and pumped to the GPC by a metering pump (model RHVOOSKY from FMI). The GPC is triggered to sample the stream every 30 min to start the analysis. Accordingly, the sample takes 30 min to travel through the chromatography columns and refractometer. At the end of the cycle, the Millenium 32 GPC software analyzes the data and the NAMW is exported as a text file for use by the control software. 3. Dynamic Modeling for the Polymerization Reactor System The reaction to be studied is the free-radical homopolymerization of MMA in a continuous stirred tank reactor (CSTR), with AIBN as the initiator and toluene as the solvent.9 The kinetic model is based on the standard free-radical polymerization mechanism shown in Figure 2 and the following standard assumptions: (i) quasi-steady-state approximation (QSSA), (ii) long-chain hypothesis (LCH), (iii) all the reaction steps are elementary and irreversible, and (iv) the rates of the reaction steps are independent of the live polymer length. This leads to the reaction rate expressions10,11

Rm ) -Cmλ0(ktm + kp) Ri ) -kiCi Rh ) (-∆Hp)kpλ0Cm Rµ0 ) (ktdλ0 + ktmCm + ktsCs)λ0 + 0.5ktcλ02 Rµ1 ) (ktmCm + ktsCs)λ1 + ktλ0λ1

(1)

where Cm, Ci, and Cs are the molar concentrations of the monomer, initiator, and solvent. The first two moments of the live polymer distributions are given by

λ1 )

(

λ0 ) (2fkiCi/kt)1/2

)

2fkiCi + (kp + ktm)λ0Cm + ktsλ0Cs Mm ktmCm + ktsCs + ktλ0

(2)

4066

Ind. Eng. Chem. Res., Vol. 41, No. 16, 2002

The rate constants appearing in eqs 1 and 2 are governed by the following set of equations:10

kl/kp ) Zl exp(-El/RT), l ) tm, ts

In the development of the reactor dynamics, it was assumed that the reactor contents are perfectly mixed, the reacting mixture has a constant heat capacity (c), and the overall heat-transfer coefficient (U) is constant. 3.2. Jacket Dynamics. The jacket dynamics are described by the following differential equation:

ktc/ktd ) Ztc exp(-Etc/RT)

dTj/dt ) RJ1(T - Tj) + RJ2(T∞ - Tj) + RJ3Q

kl ) Zl exp(-El/RT), l ) t, p, i

kt ) ktc + ktd

(3)

where T is the reactor temperature. During the conversion of monomer to polymer, there is a substantial volume contraction due to the increase in the density of the reacting mixture. For MMA, the volume contraction is more than 22% in bulk polymerization. Accordingly, this aspect must be considered when modeling the reactor. When the total volume of the reacting mixture is held constant and Fi , Fm, the flow rate of the outlet stream (F) is related to the flow rate of the inlet monomer stream (Fm) by the equation12

F ) Fm(1 + Xp)

(4)

where Q is the rate of heat addition or removal and

RJ1 )

U∞A∞ UA 1 , RJ2 ) , RJ3 ) mjcw mjcw mjcw

For these dynamics, it was assumed that the temperature of the jacket fluid inside the heating and cooling system is uniform and the heat capacity of water (cw) is constant. 3.3. Overall Dynamics. Combining the dynamic equations for the reactor and jacket gives the overall dynamic model of the reactor. The differential equations can be written in compact form as follows:

dCm/dt ) fm(Cm,Ci,µ1,T)

where Xp is the monomer conversion

dCi/dt ) fi(Ci,µ1,T) + (Ciis/V)Fi

µ1 Xp ) µ1 + MmCm

(5)

dCs/dt ) fs(Cs,µ1) + (Csis/V)Fi dµ0/dt ) fµ1(Cm,Ci,Cs,µ0,µ1,T)

is the constant volume expansion factor

(

(8)

)

CmmsMm Fm -1 ) Fm Fp

dµ1/dt ) fµ1(Cm,Ci,Cs,µ1,T) (6)

and µ1 is the mass concentration of the dead polymer chains. 3.1. Reactor Dynamics. Applying the kinetic model to a CSTR yields a set of six ordinary differential equations for the reactor dynamics. The differential equations include the mole balances on the monomer, initiator, solvent, and dead polymer (0th moment of the dead polymer distribution), a mass balance on the dead polymer (1st moment of the dead polymer distribution), and an energy balance on the reactor:

[Cmms - (1 + Xp)Cm]Fm dCm ) Rm + dt V FiCiis - (1 + Xp)CiFm dCi ) Ri + dt V dCs FmCsms + FiCsis - (1 + Xp)CsFm ) dt V dµ0 (1 + Xp)µ0Fm ) Rµ0 dt V

dT/dt ) fT(Cm,Ci,µ1,T,Tj) dTj/dt ) fTj(T,Tj) + RJ3Q

(9)

where

fm(Cm,Ci,µ1,T) ) Rm(Cm,Ci,T) + [Cmms - (1 + Xp)Cm]Fm V fi(Ci,µ1,T) ) Ri(Ci,T) fs(Cs,µ1) )

(1 + Xp)CiFm V

FmCsms - (1 + Xp)CsFm V

fµ0(Cm,Ci,Cs,µ0,µ1,T) ) Rµ0(Cm,Ci,Cs,T) -

(1 + Xp)µ0Fm V

fµ1(Cm,Ci,Cs,µ1,T) ) Rµ1(Cm,Ci,Cs,T) -

(1 + Xp)µ1Fm V

(1 + Xp)µ1Fm dµ1 ) Rµ1 dt V

fT(Cm,Ci,µ1,T,Tj) ) Rh(Cm,Ci,T)V [Tms - (1 + Xp)T]Fm + + RR(Tj - T) mc V

dT RhV [Tms - (1 + Xp)T]Fm ) + + RR(Tj - T) (7) dt mc V

fTj(T,Tj) ) RJ1(T - Tj) + RJ2(T∞ - Tj)

where

RR ) UA/mc

A list of the system parameters is given in Table 1. 3.4. Parameter Estimation. The polymerization model of eq 9 contains four unknown parameters, RR, RJ1, RJ2, and RJ3, which appear in the energy balances

Ind. Eng. Chem. Res., Vol. 41, No. 16, 2002 4067 Table 1. System Parameters Fm ) 2.78 × 10-7 m3/s V ) 1.20 × 10-3 m3 ) -0.0592 f ) 0.58 c ) 1.815 kJ/kg‚K Cmms ) 2.8046 kmol/m3 Csms ) 6.5715 kmol/m3 Ciis ) 0.2588 kmol/m3 Csis ) 9.3878 kmol/m3 Tms ) 296.6 K T∞ ) 296.6 K Mm ) 100.12 kg/kmol Fm ) 936 kg/m3 Fp ) 1166 kg/m3

Table 2. Estimated Parameters RR ) 0.005 655 s-1 RJ1 ) 0.000 559 7 s-1

Zt ) 9.800 × 107 m3/(kmol‚s) Et ) 2.933 × 103 kJ/kmol Zp ) 4.917 × 105 m3/(kmol‚s) Ep ) 1.820 × 104 kJ/kmol Zi ) 1.053 × 1015 s-1 Ei ) 1.284 × 105 kJ/kmol Ztm ) 4.661 × 109 m3/(kmol‚s) Etm ) 7.418 × 104 kJ/kmol Zts ) 1.010 × 103 m3/(kmol‚s) Ets ) 4.769 × 104 kJ/kmol Ztc ) 3.956 × 10-4 m3/(kmol‚s) Etc ) -1.711 × 104 kJ/kmol -∆Hp ) 57 800 kJ/kmol

of the reactor and jacket. Therefore, these parameters must be estimated to complete the model and construct the multirate controller. The estimation procedure for each parameter involves standard least-squares analysis (linear regression) of experimental data fit to equations derived from the original continuous-time model: (i) Estimation of rR. Under the condition of only solvent in the reactor, the energy balance of the reactor becomes dT/dt ) RR(Tj - T). When a step change in the jacket temperature occurs, the dynamics of the reactor temperature follow the equation

t)-

(

T - Tcw 1 ln RR T(0) - Tcw

)

(10)

d(Tj - T∞)/dt ) RJ3P - RJ3(UA + U∞A∞)(Tj - T∞) (11) When a step change in the heater power (P) is applied, the following equation results:

-

(

(iv) Calculation of rJ1 and rJ2. The straightforward calculation of RJ1 and RJ2 is as follows: RJ1 ) UARJ3 ) 0.000 559 4 s-1 and RJ2 ) U∞A∞RJ3 ) 0.000 452 1 s-1. Table 2 summarizes the results for the estimated parameters. 3.5. Discrete-Time Model. The digital implementation of the nonlinear multirate model-algorithmic controller requires a discrete-time process model. The continuous-time model of eq 9 was discretized using Euler’s method to give the following overall discretetime model:

Cm(k+1) ) Cm(k) + fm[Cm(k),Ci(k),µ1(k),T(k)]∆t Ci(k+1) ) Ci(k) + fi[Ci(k),µ1(k),T(k)]∆t + Cs(k+1) ) Cs(k) + fs[Cs(k),µ1(k)]∆t +

Ciis∆t Fi(k) V

Csis∆t Fi(k) V

µ0(k+1) ) µ0(k) + fµ0[Cm(k),Ci(k),Cs(k),µ0(k),µ1(k),T(k)]∆t µ1(k+1) ) µ1(k) + fµ1[Cm(k),Ci(k),Cs(k),µ1(k),T(k)]∆t

The step change in the jacket temperature is achieved by changing the cooling water valve from the closed position to the fully open position with no power to the heater. The slope of the fitted regression line to ln[(T Tcw)/(T(0) - Tcw)] data for different time instants gives the estimate of RR ) 0.005 655 s-1, which yields UA ) 0.010 36 kJ/(K‚s). (ii) Estimation of U∞A∞. Under the conditions of only solvent in the reactor, steady-state operation (T and Tj constant), and no inlet cooling water (Fcw ) 0), the dependence of the difference between the jacket temperature and the ambient temperature (Tj - T∞) on the heater power (P) is linear. In this case, the slope of the fitted regression line to Tj - T∞ data for different settings of P gives the estimate of U∞A∞ ) 0.008 368 kJ/(K‚s). (iii) Estimation of rJ3. Under the conditions of only ambient air in the reactor and no inlet cooling water (Fcw ) 0), the energy balance for the jacket becomes

t)

RJ2 ) 0.000 452 1 s-1 RJ3 ) 0.054 03 K/kJ

)

Tj - T∞ - P/(UA + U∞A∞) 1 ln RJ3(UA + U∞A∞) Tj(0) - T∞ - P/(UA + U∞A∞) (12)

The step change is achieved by changing the heater power from 0 to 1.5 kJ/s. Using the above equation, the slope of the fitted regression line to ln{[Tj - T∞ - P/(UA + U∞A∞)]/[Tj - T∞ - P/(UA + U∞A∞)]} data gives the estimate of RJ3 ) 0.054 03 K/kJ.

T(k+1) ) T(k) + fT[Cm(k),Ci(k),µ1(k),T(k),Tj(k)]∆t Tj(k+1) ) Tj(k) + fTj[T(k),Tj(k)]∆t + RJ3∆tQ(k) (13) The time step ∆t was chosen to be 30 s (0.5 min). The manipulated inputs of the system are the inlet initiator steam flow rate (u1 ) Fi) and the net rate of heat addition or removal from the jacket (u2 ) Q). The input actuation period is equal to the discretization time step ∆t. The controlled outputs of the system are the NAMW (y1 ) µ1/µ0) and the reactor temperature (y2 ) T). The NAMW is measured every 30 min with a 30 min time delay, while the reactor temperature is measured every 0.5 min without delay (φ2 ) 0). Because the calculated input moves are actuated every 0.5 min, N1 ) φ1 ) 60 and N2 ) 1. 3.6. Dynamic Characteristics of the Polymerization Reactor System. The dynamic model developed in the previous subsections is a two-input seven-state two-output discrete-time system with measurement delays. Before the multirate MAC algorithm of part 18 can be applied, it is important to verify that the requirements for applicability of the theory are indeed satisfied. First, the requirement of finite relative orders must be checked. A straightforward calculation shows that the relative orders of the outputs are r1 ) r2 ) 2, excluding the measurement delays. Second, reference steady-state conditions must be chosen for the experiments, around which the dynamic characteristics of the process must be evaluated. Using the steady-state equations from the model, it was found that, for Fi ) 1.355 × 10-8 m3/s and Q ) 0.399 kJ/s, the corresponding steady-state values for the reactor temperature and NAMW are T ) 343 K and NAMW ) 45 000 kg/kmol. These represent reasonable set-point values for the

4068


Table 3. Reference Steady-State Conditions Cm ) 2.0583 kmol/m3 Ci ) 0.011 36 kmol/m3 Cs ) 7.1468 kmol/m3 T ) 343.0 K

Tj ) 343.6096 K µ0 ) 0.001 769 kmol/m3 µ1 ) 79.5888 kg/m3

Table 4. Eigenvalues of JOL and JINV eigenvalues of JOL eigenvalues of JINV

0.993 165 0.989 358 ( 0.001 542 4i 0.993 159

0.819 764 0.993 455 0.987 405

0.993 164 0.992 289 0.993 436

Table 5. Parameters for Coordination Rules cwFw ) 4200 kJ/(m3‚K) Fcwmax ) 2.51 × 10-4 m3/s Tcw ) 290.4 K

Pmax ) 3.0 kJ/s Fimax ) 4.0 × 10-8 m3/s

system outputs. Table 3 provides the complete set of values for the state variables corresponding to this steady state. Using these steady-state values, the Jacobians JOL (for the open-loop system) and JINV (for the inverse of the delay-free part) have been calculated along with their eigenvalues. The results are summarized in Table 4. Because all eigenvalues are inside the unit disk, the system is both hyperbolically stable and hyperbolically minimum phase around the reference steady state. Consequently, the application of the multirate MAC controller developed in part 1 will lead to an asymptotically stable closed-loop system as long as the reference trajectories are stable (see theorem 3 of part 1). Finally, the speed of the process response can be estimated by simulating the output responses to step changes in the inputs. Parts a and b of Figure 3 show the open-loop output responses for step changes in the initiator flow rate Fi and the rate of heat addition to the reactor jacket Q, respectively. For Figure 3a, the initiator flow rate was changed to 1.153 × 10-8 m3/s from its reference steady-state value, while for Figure 3b, the rate of heat addition was changed to 0.379 kJ/s. It is observed that, for either step change in the input, it takes a few hours for each output to settle to its new steady state. Thus, the two outputs appear to have a roughly equal speed of response. 3.7. Practical Considerations. As with many physical systems, active input constraints are present in the operation of the reactor. Constraints appear in both the initiator flow rate (u1 ) Fi) and the net rate of heat addition or removal from the reactor jacket (u2 ) Q). Because there are two actual manipulated inputs driving the heating and cooling of the jacket (P and Fcw), they must be coordinated to combine both of the inputs into a single manipulated input such that Q ) f(P,Fcw). One possible approach is to combine the manipulated inputs linearly. The disadvantage associated with a linear combination is that the heating and cooling almost always take place simultaneously, making the operation undesirable energywise. Therefore, the following coordination approach is used to combine the two actual manipulated inputs.14 First, combine the actual manipulated inputs such that

u2 ) P - FcwcwFw(Tj - Tcw)

(14)

where cw, Fw, and Tcw are the heat capacity, density, and temperature of the cooling water, respectively. After u2 has been calculated from the control algorithm,

Figure 3. Output responses for (a) a step change in Fi and (b) a step change in Q.

P and Fcw are set according to the following coordination rules:

{

{

u2 if 0 e u2 < Pmax P P) max if u2 g Pmax if u2 < 0 0

Fcw )

-u2

(15)

if -FcwmaxcwFw(Tj - Tcw) e u2 < 0

cwFw(Tj - Tcw) Fcwmax

if u2 < -FcwmaxcwFw(Tj - Tcw)

0

if u2 g 0 (16)

where Pmax and Fcwmax are the maximum heater power and coolant flow rate. In addition, the initiator flow rate (u1 ) Fi) is bounded by the pump capacity. The actual manipulated input Fi is calculated from

{

0

if u1 < 0

Fi ) u1 if 0 e u1 < Fimax Fimax if u1 g Fimax

(17)

where u1 is calculated from the multirate controller and Fimax is the maximum initiator flow rate. It is noted that the actual values of the manipulated inputs given above are used to drive the state estimator equations and not the calculated values from the control algorithm. The parameters used in the coordination rules are given in Table 5. 3.8. Controller Implementation. Because the reactor is a two-input two-output system, the multirate MAC controller developed in part 18 will involve two nonlinear algebraic equations that must be solved online to calculate the control action. For relative orders r1 )

Ind. Eng. Chem. Res., Vol. 41, No. 16, 2002 4069 Table 6. Initial Conditions for Reactor Startup Cm(0) ) 2.8046 kmol/m3 Ci(0) ) 0 kmol/m3 Cs(0) ) 6.5715 kmol/m3 T(0) ) 296.6 K

Tj(0) ) 296.6 K µ0(0) ) 0 kmol/m3 µ1(0) ) 0 kg/m3

r2 ) 2, N2 ) 1, and φ2 ) 0, the following equations result from eq 56 of part 1: 1

[

([ ] )

h1 [Φ[xM(k),u(k)]] ) (1 - R1) y1sp - y1

[ ([ ]

h1 xM

)]]

k N - φ1 N1 1

k N + N1 1

+ R1h11[xM(k)]

h21[Φ[xM(k),u(k)]] ) (1 - R2)[y2sp - y2(k) + h2[xM(k)]] + R2h21[xM(k)] (18) where xM is the vector of model states xM(k+1) ) Φ[xM(k),u(k)], which is obtained by online simulation of the process model (13). When the necessary function compositions are performed in eq 18, it is found that the first equation is only a function of the input u1 and independent of u2. Therefore, the two algebraic equations can be solved separately. The standard secant method for finding roots of nonlinear equations was used to solve for u1. The second equation is linear in u2, and therefore the solution for u2 is readily obtainable once u1 is found. Finally, the coordination rules (15)(17) are used to compute the actual manipulated inputs Fi, P, and Fcw from u1 and u2. 4. Experimental Procedure After monomer purification by standard methods,13,14 the reactor is loaded with 1.2 × 10-3 m3 of the monomer solution (30% methyl methacrylate-70% toluene by volume at room temperature). The reactor feed bottles are then loaded with the appropriate solutions. The monomer feed solution has the same composition as the reactor solution, and the initiator feed solution has a concentration of 0.2588 (kmol AIBN)/m3. Aldrich Chemical Co. supplied all HPLC-grade toluene, MMA, and AIBN used in the experiments. When the loading of the solutions is complete, the system is purged of oxygen (a reaction inhibitor) by bubbling nitrogen through it for 1 h prior to the start of the reaction. A nitrogen environment is also maintained over the solutions in the reactor and feed bottles during the course of an experiment. To start an experimental run, the computer program is first initialized to perform the 2 h reactor startup. The initial conditions are given in Table 6. During the startup, the initiator flow rate is set to its maximum value of 4 × 10-8 m3/s to bring the NAMW down to its set-point value in the shortest possible time. When the NAMW reaches a value within 10 000 kg/kmol of the set point, the values calculated from the controller equations are then used to adjust the initiator flow rate. It is noted that the GPC does not take its first sample until an hour after startup to allow the NAMW to reach values within the molecular weight resolution of the GPC columns. The reactor temperature is brought up to the operating temperature according to the values of the heater power calculated from the controller equations. At each sampling period (0.5 min), the computer program checks for new NAMW data that were output from the GPC software in addition to receiving the

reactor temperature. If there are new NAMW data, they are incorporated into the controller calculations until the next data point is received. In the event of a fault in the analysis of the MWD by the GPC software, the computer program uses the last available measurement in the controller calculations. 5. Controller Performance The performance of the multirate model-algorithmic controller is investigated in terms of steady-state regulation and tracking of step changes in the set points. Through the open-loop simulations of the polymerization model (see section 3.6), it was found that the speed of response for both outputs is approximately equal when a step change in a manipulated input is applied. The approximate time constant for a step change in the heater power (P) is slightly less than 1 h for both outputs, whereas the approximate time constant for a step change in the initiator flow rate (Fi) is about 2 h. This observation suggests that it is reasonable to request equally fast closed-loop responses for both outputs. Therefore, the same values of R will be used for both outputs (R ) R1 ) R2). This leaves the control law with only one adjustable parameter, which makes the tuning of the controller significantly easier. The experimental results begin with the steady-state operation of the reactor, followed by step changes in the NAMW and reactor temperature set points. 5.1. Steady-State Operation. Parts a and b of Figure 4 show the NAMW and reactor temperature profiles, respectively, for constant set-point values. For this experimental run, the NAMW set point is 45 000 kg/kmol and the temperature set point is 343.0 K, with R1 ) R2 ) 0.95. As is seen in Figure 4a, the controller capably regulates the NAMW to well within 500 kg/kmol of the requested set point. Similarly, Figure 4b shows that the reactor temperature is sufficiently controlled to 0.5 K of the set-point value. The corresponding initiator flow rate (Fi) profile is shown in Figure 4c, and Figure 4d shows the heater power (P) and coolant flow rate (Fcw) profiles. The effect of tight input constraints in the presence of infrequent measurements is clearly seen in the initiator flow-rate profile. When a new measurement is received, the controller aggressively compensates for the error between the model and measurement. If the NAMW measurement is higher than the model value, the initiator flow rate jumps to its upper constraint before going toward the steady-state value. The addition of initiator will increase the rate of chain initiation reactions, which creates new shorter chains and effectively decreases the NAMW. Analogously, a lower NAMW measurement causes the initiator flow rate to jump to its lower constraint. This decreases the rate of chain initiation reactions, which increases the NAMW. It is also seen that the heater power P varies slightly around its steady-state value. Because u2 is positive, the coolant flow rate is zero (Fcw ) 0) for the constant set-point experiment. 5.2. Step Change in the NAMW Set Point. Step changes in the NAMW set point are now considered with varying controller parameters, R (R1 ) R2 for all cases). Figures 5-7 show the profiles for a step change in the NAMW set point from 45 000 to 50 000 kg/kmol with a constant temperature set point of 343.0 K. For traditional industrial controllers, a change in the molecular weight of a polymer (grade change) might require the shutdown and restart of the reactor5 or at

4070


Figure 4. (a) NAMW profile. (b) Temperature profile. (c) Initiator flow-rate profile. (d) Heater power and coolant flow-rate profiles for constant set points, R ) 0.95.

Figure 5. (a) NAMW profile. (b) Temperature profile. (c) Initiator flow-rate profile. (d) Heater power and coolant flow-rate profiles for a step change in the NAMW set point, R ) 0.975.

least the retuning of the controllers. Figures 5a, 6a, and 7a show that, with the nonlinear multirate controller, a step change in the molecular weight is easily achieved regardless of controller parameters and without the need for controller retuning. It is seen that decreasing the R value does not significantly increase the speed of the NAMW response. This is due to the lower constraint on the initiator flow rate. Figures 5c, 6c, and 7c show that the initiator flow rate initially drops to zero. As previously mentioned, a decrease in the initiator flow rate will increase the NAMW. A controller with a smaller R will calculate a larger negative flow rate of

the initiator. Because it is physically impossible to remove the initiator from the reactor at any rate, the speed of the NAMW response is limited. However, the increase in the aggressiveness of the controller with a smaller R is clearly seen in the manipulated input profiles. As R decreases, the initiator flow rate spends more time at the constraints when a measurement is received in an effort to eliminate the error in the shortest time possible. It is also seen in Figures 5d, 6d, and 7d that a decrease in R leads to increased variations in the heater power around its operating value. Again, the controller is making more aggressive input moves


Figure 6. (a) NAMW profile. (b) Temperature profile. (c) Initiator flow rate profile. (d) Heater power and coolant flow-rate profiles for a step change in the NAMW set point, R ) 0.95.

Figure 7. (a) NAMW profile. (b) Temperature profile. (c) Initiator flow-rate profile. (d) Heater power and coolant flow-rate profiles for a step change in the NAMW set point, R ) 0.90.

to eliminate the error. The effects are seen in Figures 5b, 6b, and 7b, where a smaller R decreases the drift of the reactor temperature. The experiments for a NAMW set-point change can be summarized as follows: (i) Decreasing the controller parameter R does not significantly increase the speed of the NAMW response due to the constraints on the initiator flow rate. (ii) Decreasing the controller parameter R improves control of the reactor temperature during a set-point change in the NAMW.

(iii) Decreasing the controller parameter R causes the initiator flow rate to remain longer at the constraints when a slowly sampled measurement is received. 5.3. Step Change in the Temperature Set Point. Figures 8-10 show the profiles for a step change in the reactor temperature set point from 343.0 to 348.0 K, with a constant NAMW set point of 45 000 kg/kmol. It is seen in Figures 8b, 9b, and 10b that the temperature response becomes faster with decreasing R. It takes 1 h for the temperature to reach the new set point with R

4072


Figure 8. (a) NAMW profile. (b) Temperature profile. (c) Initiator flow-rate profile. (d) Heater power and coolant flow-rate profiles for a step change in the temperature set point, R ) 0.975.


) 0.975 but only 0.3 h with an R ) 0.90. The corresponding heater power profiles are shown in Figures 8d, 9d, and 10d. Initially, the heater power jumps to 0.86, 1.30, and 2.21 kJ/s for R ) 0.975, 0.95, and 0.90, respectively, before proceeding to the new operating value. The initial heater power increases with decreasing R to force the faster temperature response. The corresponding NAMW profiles are seen in Figures 8a, 9a, and 10a. As the temperature response becomes

faster, the error in the NAMW becomes larger. An increase in temperature will increase the rate of chain termination reactions, decreasing the NAMW. As shown in Figures 8c, 9c, and 10c, the controller initially compensates for the increase in temperature by stopping the flow of initiator into the reactor. As the change in temperature becomes faster, the controller cannot fully compensate for the temperature effect due to the lower initiator constraint. For R ) 0.975 the controller can



keep the NAMW at the set point, while for R ) 0.90 the NAMW drops to nearly 43 000 kg/kmol and takes 2 h to recover. The experiments for a temperature set-point change can be summarized as follows: (i) Decreasing the controller parameter R increases the speed of the reactor temperature response. (ii) Decreasing the controller parameter R deteriorates the control of the NAMW during a set-point change in the reactor temperature due to the constraints on the initiator flow rate. It is observed that for a step change in the NAMW set point the outputs are decoupled for all values of R, but for a step change in the temperature set point the outputs are only decoupled for larger values of R. The inadequate decoupling is due to the lower constraint on the initiator flow rate. This poses a design tradeoff question of whether the speed of response or output decoupling is more important in the process operation. If the product specifications have a broad tolerance for the molecular weight of the polymer, then it would be advantageous from a disturbance rejection perspective to speed up the controller with a smaller R. If a narrow tolerance is sought for the molecular weight, then a larger R should be used because of its decoupling effect from the reactor temperature. 6. Conclusion This paper experimentally implements the nonlinear multirate model-algorithmic controller8 to regulate the reactor temperature and NAMW of a continuous-stirred tank polymerization reactor. A discrete-time model of the free-radical MMA polymerization was presented and used for the construction of the multirate control algorithm. By variation of the tunable controller parameters, the performance of the controller was studied in the presence of active input constraints. It was shown that the controller adequately tracks step changes in the output set points. Decreasing the tunable controller

parameters deteriorated the control of the NAMW output but improved the control of the temperature output. This deteriorated control is due to the lower input constraint on the inlet initiator flow rate. Overall, the proposed controller represents an efficient and effective method to deal with the difficult problem of multirate output measurements and to achieve smooth grade transitions. Acknowledgment Financial support from the National Science Foundation through Grant CTS-9403432 is gratefully acknowledged. Notation A, A∞ ) reactor jacket and ambient jacket heat-transfer surface areas, m2 c ) heat capacity of the reacting mixture, kJ/(kg‚K) Ci ) concentration of the initiator, kmol/m3 Ciis ) concentration of the initiator in the inlet initiator stream, kmol/m3 Cm ) concentration of the monomer, kmol/m3 Cmms ) concentration of the monomer in the inlet monomer stream, kmol/m3 Cs ) concentration of the solvent, kmol/m3 Csis ) concentration of the solvent in the inlet initiator stream, kmol/m3 Csms ) concentration of the solvent in the inlet monomer stream, kmol/m3 cw ) heat capacity of water, kJ/(kg‚K) Ei ) activation energy for ki, kJ/kmol Et, Ep ) activation energies for kt and kp, kJ/kmol Etc ) activation energy for ktc, kJ/kmol Etm, Ets ) activation energies for ktm and kts, kJ/kmol f ) initiator efficiency Fcw ) flow rate of inlet cooling water, m3/s Fcwmax ) maximum flow rate of inlet cooling water, m3/s Fi ) initiator flow rate, m3/s

4074


Fimax ) maximum initiator flow rate, m3/s Fm ) monomer flow rate, m3/s ki ) rate constant for initiation reactions, s-1 kp, kt ) rate constants for propagation and overall termination reactions, m3/(kmol‚s) ktc, ktd ) rate constants for termination reactions by combination and disproportionation, m3/(kmol‚s) ktm, kts ) rate constants for chain-transfer reactions to monomer and solvent, m3/(kmol‚s) m ) mass of the reacting mixture inside the reactor, kg mj ) mass of fluid in the heating/cooling system, kg Mm ) molecular weight of the monomer, kg/kmol Ni ) ratio of the output sampling period to the input actuation period P ) heater power, kJ/s Pmax ) maximum heater power, kJ/s Q ) net rate of heat addition or removal, kJ/s R ) universal gas constant, 8.314 kJ/(kmol‚K) ri ) relative order of the output T ) reactor temperature, K T∞ ) ambient temperature, K Tcw) cooling water temperature, K Tj ) reactor jacket temperature, K Tms ) inlet monomer stream temperature, K u ) manipulated input vector U, U∞ ) reactor jacket and ambient jacket overall heattransfer coefficients, kJ/(s‚m2‚K) V ) volume of the reacting mixture, m3 x ) process state vector Xp ) conversion y ) process output vector Zi ) frequency factor for ki, s-1 Zt, Zp ) frequency factors for kt and kp, m3/(kmol‚s) Ztc ) frequency factor for ktc, m3/(kmol‚s) Ztm, Zts ) frequency factors for ktm and kts, m3/(kmol‚s) Greek Letters R ) tunable controller parameter RJ1, RJ2 ) heat-transfer constants for the reactor jacket, s-1 RJ3 ) heat-transfer constant for the reactor jacket, K/kJ RR ) heat-transfer constant for the reactor, s-1 -∆Hp ) heat of propagation reactions, kJ/kmol ∆t ) sampling period, s ) volume expansion factor φ ) measurement time delay λ0, λ1 ) 0th and 1st moments of the radical distributions, kmol/m3 and kg/m3 µ0, µ1 ) 0th and 1st moments of the dead polymer distributions, kmol/m3 and kg/m3

Fm, Fp ) densities of monomer and polymer, kg/m3 Fw ) density of water, kg/m3

Literature Cited (1) Nunes, R. W.; Martin, J. R.; Johnson, J. F. Influence of Molecular Weight and Molecular Weight Distribution on Mechanical Properties of Polymers. Polym. Eng. & Sci. 1982, 4, 205. (2) 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. (3) Adebekun, D. K.; Schork, F. J. Continuous Solution Polymerization Reactor Control. 2. Estimation and Nonlinear Reference Control during Methyl Methacrylate Polymerization. Ind. Eng. Chem. Res. 1989, 28, 1846. (4) Mutha, R. K.; Cluett, W. R.; Penlidis, A. On-line Nonlinear Model-Based Estimation and Control of a Polymer Reactor. AIChE J. 1997, 43, 3042. (5) Ogunnaike, B. A. On-Line Modeling and Predictive Control of an Industrial Terapolymerization Reactor. Int. J. Control 1994, 59, 711. (6) Ohshima, M.; Hashimoto, I.; Ohno, H.; Takeda, M.; Yoneyama, T.; Gotoh, F. Multirate Multivariable Model Predictive Control and its Application to a Polymerization Reactor. Int. J. Control 1994, 59, 731. (7) Tatiraju, S.; Soroush, M.; Ogunnaike, B. A. Multirate Nonlinear State Estimation with Application to a Polymerization Reactor. AIChE J. 1999, 45, 769. (8) Niemiec, M.; Kravaris, C. Nonlinear Model-Algorithmic Control. 1. Theory. Ind. Eng. Chem. Res. 2002, 41, 4054-4063. (9) Soroush, M.; Kravaris, C. Multivariable Nonlinear Control of a Continuous Polymerization Reactor: An Experimental Study. AIChE J. 1993, 39, 1920. (10) Baillagou, P. E.; Soong, D. S. Molecular Weight Distribution of Products of Free Radical Nonisothermal Polymerization with Gel Effect. Simulation for Poly(methyl methacrylate). Chem. Eng. Sci. 1985, 40, 87. (11) Ray, W. H. On the Mathematical Modeling of Polymerization Reactors. J. Macromol. Sci., Rev. Macromol. Chem. 1972, C8, 1. (12) Schmidt, A. D.; Ray, W. H. The Dynamic Behavior of Continuous Polymerization Reactors: I. Isothermal Solution Polymerization in a CSTR. Chem. Eng. Sci. 1981, 36, 1401. (13) Collins, E. A.; Bares, J.; Billmeyer, F. W. Experiments in Polymer Chemistry; Wiley: New York, 1973. (14) Soroush, M.; Kravaris, C. Nonlinear Control of a Batch Polymerization Reactor: An Experimental Study. AIChE J. 1992, 38, 1429.

Received for review June 28, 2001 Revised manuscript received January 7, 2002 Accepted April 5, 2002 IE010560M

Nonlinear Multirate Model-Algorithmic Control. 2. Experimental

Recommend Documents