Multivariable Control of a Simulated Industrial Gas-Phase

Ind. Eng. Chem. Res. 2003, 42, 2349-2364

2349

Multivariable Control of a Simulated Industrial Gas-Phase Polyethylene Reactor E. Ali,* K. Al-Humaizi, and A. Ajbar Chemical Engineering Department, King Saud University, P.O. Box 800, Riyadh 11421, Saudi Arabia

A multivariable control problem of an industrial gas-phase polyethylene reactor is investigated in this paper. This multivariable control problem affects directly the polymer properties and consequently the plant economics. The control task is particularly challenging because of the two-time-scale behavior of the process and because of the multirate sampling of the process controlled variables. Two control schemes and two control algorithms are tested and compared. One control scheme considers the overall control problem using one multivariable control algorithm such as a model predictive controller (MPC). Another control scheme separates the overall control problem into two levels. One level takes care of the controlled variables with fast dynamics via a multiloop proportional-integral algorithm, and the other level takes care of the slow controlled variables via a multivariable control algorithm such as MPC. The simulated comparison revealed the superiority of the second scheme only in providing a faster response. Other than the speed of response, the performances of the two schemes are almost comparable. Further, the first scheme may outperform the second one in the presence of input limitation. On the other hand, linear MPC and nonlinear MPC algorithms were compared for the two control schemes. Simulation results indicated that the nonlinear MPC outweighs the linear MPC in terms of performance, tuning requirements, and robustness. Introduction Polyethylene is considered to be the world largest produced synthetic commodity polymer. During the past decade, consistency requirements in the polymer properties have led to a growing interest in product-quality control schemes for polymerization reactors. Control of polymerization reactors has long been known to be a difficult task because of the high nonlinearity of the reactions involved and the strong interaction between the reactor variables. However, there is relatively little work on the control of gas-phase polymerization of ethylene in fluidized-bed reactors. Industrial gas-phase fluidized-bed polyethylene reactors are operated in a narrow temperature range between 75 and 130 °C, and without appropriate temperature control, they are known to be prone to unstable behavior, temperature oscillations, and incursions toward serious runaway.1-4 Stabilization of polyethylene reactors is therefore a challenging problem and needs to be addressed through good control. Dadebo et al.5 and Ali et al.6,7 have addressed the stabilization of the reactor temperature through different approaches. Recently, Seki et al.8 have also studied the stabilization of gas-phase polyethylene reactors through an adequate tuning of the proportionalintegral-derivative controller. All of the aforementioned studies dealt only with the tight regulation of the bed temperature to ensure reactor stability. McAuley and MacGregor9 have studied the control of the polymer quality through manipulating the feed flows. In this paper, we deal with a more comprehensive control problem. Specifically, the control of the reactor temperature and pressure in addition to the gas partial pressures is considered in this paper. This control objective is more appealing because maintaining * To whom correspondence should be addressed. Tel.: +9961-467-6871. Fax: +966-1-467-8770. E-mail: amkamal@ksu. edu.sa.

the ratio of the partial pressures at specific values has a significant role in the polymer quality that directly affects plant economics. This multivariable control problem is, however, challenging because strong crossloop interaction prevails. To add to the difficulty, the process is well-known to have complex dynamics due to the feedback interaction induced by the recycle stream.10,11 In general, the overall process with recycle may exhibit a time-scale separation in its dynamics. Specifically, the dynamics of the recycle network can be segregated into fast-time and slow-time control loops.12 This phenomenon of two time scales creates a significant detrimental effect on the overall multivariable control system. The differential time scale is paramount in many other chemical processes where the ratio of the residence times of individual units is very large. The situation may lead to controller ill-conditioning, stiffness of the numerical simulation, and closedloop instability.13 Another problem that is related to the multivariable control of the polyethylene reactor is the multirate sampling because some controlled variables are measured less frequently than others. This paper deals with multivariable control of the polyethylene reactor using two approaches. The first approach is to consider the control problem as one multivariable control system. In this case, a multivariable control algorithm such as model predictive control (MPC)14 is used. The second approach12 is to consider two separate control levels. One level is a multivariable system consisting of the controlled variables with slow dynamics, while the other level is a multiloop singleinput single-output (SISO) system containing the controlled variables with fast dynamics. For the first level, an MPC algorithm is used, while a standard proportional-integral (PI) algorithm is used for the second level. A unified sampling time will be used for all loops involved in the multivariable MPC system. The sampling time can be set as small as that used for the slowest measured variable. For the SISO PI loop, a different sampling time can be used, which can be set

10.1021/ie020386h CCC: $25.00 © 2003 American Chemical Society Published on Web 04/23/2003

2350

Ind. Eng. Chem. Res., Vol. 42, No. 11, 2003

and is given below. This model is chosen because its kinetic parameters were validated against plant data.15 The definitions of the various states and parameters of the model are given in the Nomenclature section.

Vg

dCM1 ) FM1 - xM1Bt - RM1 dt

(1)

Vg

dCM2 ) FM2 - xM2Bt - RM2 dt

(2)

dCM2 ) FH - xHBt - RH dt

(3)

dCN ) FN - xNBt dt

(4)

Vg

Vg

dYc ) Fcac - kdYc - OpYc/Bw dt (MrCp,r + BwCp,p)

Figure 1. Schematic of the polyethylene reactor.

as small as that used for the fastest measured variable. The limitation on the fastest sampling rate that can be used and the conflict caused by the variable dynamics response of the different loops make the control task more challenging. Therefore, one objective of the paper is to study the relative merits of both approaches with respect to the two-time-scale behavior. Another objective of the paper is to compare the effectiveness of two multivariable control algorithms, namely, linear MPC (LMPC) based on the finite impulse response (FIR) model and nonlinear MPC (NLMPC). In most industrial applications, LMPC is more preferable because of its modest requirement for modeling and computations. For example, unlike the nonlinear theoretical models used by NLMPC, FIR models can be directly extracted from plant input/output data. Moreover, LMPC law is generated from the solution of, at the most, a quadratic programming problem, while the NLMPC requires the solution of a nonlinear optimization problem, which is computationally demanding and usually slow for real-time control. For these reasons, the added features of NLMPC are not yet justified. However, the robustness issue is more paramount in LMPC than in NLMPC because the average FIR model does not provide the required performance over a wide range of operating conditions. In fact, the LMPC performance may substantially degrade for open-loop unstable or slow dynamics processes. For this reason, we will present a disturbance estimation method to improve the performance of LMPC in these cases. Furthermore, the LMPC performance is severely affected by the ill-conditioning of the controller gain matrix. Careful tuning of the LMPC is necessary to obtain stable responses. Therefore, we present an ad hoc method to reduce illconditioning and to improve the LMPC performance. Process Model The polyethylene reactor process is depicted in Figure 1. The process model was developed by McAuley et al.2

MgCp,g

(5)

dT ) HF + HG - HR - HT - HP dt (6)

dTg ) FgCp,g(Tgi - Tg) + FwCp,w(Twi - Two) dt (7) Pt ) (CM1 + CM2 + CH + CN)RT Tgi )

(

Pt Pt + ∆P

)

(8)

(1-η)/η

T

(9)

FwCp,w(Twi - Two) ) 0.5UA[(Two + Twi) - (Tgi + Tg)] (10) where

HF ) (FM1Cp,M1 + FM2Cp,M2 + FHCp,H + FNCp,N) (Tf - Tref) (11) HG ) FgCp,g(Tg - Tref)

(12)

HT ) (Fg + Bt)Cp,g(T - Tref)

(13)

HP ) OpCp,p(T - Tref)

(14)

HR ) Mw1RM1∆Hr

(15)

Op ) Mw1RM1 + Mw2RM2

(16)

RM1 ) CM1Yckp1e-(E/R)(1/T-1/Tref)

(17)

RM2 ) CM2Yckp2e-(E/R)(1/T-1/Tref)

(18)

Cp,g )

∑xiCp,i

(19)

The model equations listed above are slightly modified from those given by McAuley et al.2 For simplicity, only one active site for the catalyst is considered here. Moreover, the energy balance around the cooler considers the dynamics of the recycle temperature explicitly instead of the heat removal as used by McAuley et al.2 In due course, the cooling process is modeled as a wellmixed system. The thermal effect of the recycle com-

Ind. Eng. Chem. Res., Vol. 42, No. 11, 2003 2351 Table 1. Steady-State Operating Condition CM1 CM2 CH CN Yc T Bt Fw

297.06 mol/m3 116.17 mol/m3 105.78 mol/m3 166.23 mol/m3 5.849 mol 82.7 °C 10.39 mol/s 5.6 × 105 g/s

FM1 FM2 FH FN Tg Two Twi Fc

131.13 mol/s 3.5100 mol/s 1.6000 mol/s 2.5200 mol/s 324.7 K 308 K 293 K 2 kg/s

pressor is also included in this model. Note that the partial pressure of the reactants can be calculated directly form the reactant concentrations using an ideal gas law. In eq 9, η is the isentropic factor and is set to 0.5 in this paper to make the model predictions match the industrial data. The modified model was compared to the original model to ensure that its static and dynamic characteristics are preserved. It should be noted that the reactor model (eqs 1-19) will be used to simulate the plant in the various control tests shown in this paper. The same model will also be used in the NLMPC algorithm to compute the control actions except some changes in the values of the model parameters as will be described later. On the other hand, a linearization of the reactor model in the form of FIR will be used in the LMPC algorithm. Plant Operating Condition. In addition to the above modification, the utilization of the model is more comprehensive in this work. Although McAuley et al.2 have determined all of the process parameters, they have not determined the steady-state operation condition for the entire process variables. The steady-state values for CM2, CH, CN, Yc, Two, Tg, FM2, FN, and FH were not determined. This is because the authors focused their analysis only on the reactor temperature and the monomer concentration. In this work, the objective is to control the partial pressures of the gases, the total pressure, and the reactor temperature. Therefore, an initial steady state for the entire process should be determined. One operating condition can be determined from solving the following steady-state optimization problem:

min

CM1,CM2,CH,CN,Yc,T,FM2,FH,FN,Bt,Two,Tg

remaining variables and parameters that are not included in the set of design parameters in the above minimization problem are fixed at their nominal values given by McAuley et al.2 The results of the optimization problem are listed in Table 1, while the values of the other fixed parameters are given in Table 2. To be able to make use of the dynamic model, the value of the holdup inside the cooler, Mg, needs to be determined. The value is assumed to be around 2 × 106 mol, which is estimated roughly from the cooler residence time given by McAuley et al.2 and from the coolant flow rate. Online NLMPC Algorithm In this paper, the structure of the MPC version developed by Ali and Zafiriou14 that utilizes directly the nonlinear model for output prediction is used. The proposed NLMPC scheme belongs to the nonlinear programming (NLP) based algorithms. Sequential quadratic programming is used to solve these NLP problems. The solution is carried out sequentially in two loops. In the outer loop, an optimization routine iteratively selects the new values for the manipulated variables (MVs). In the inner loop, an ODE solver is used to integrate the model equations at each iteration of the optimization. A usual MPC formulation solves the following online optimization: P

min

M

||Λ∆u(tk+i-1)||2 ∑ i)1

||Pt -

AT∆U(tk) e b

f(X h ,U h) ) 0

(21)

The process operating conditions, i.e., the values of the parameters shown in the objective function, are determined such that the monomer partial pressure and the total pressure are fixed at the desired industrial values (P /M1 and P /t ), while satisfying the algebraic constraints, f(X h ,U h ), that represent the steady-state model. The desired values for PM1 fall in the range of 7-9 atm, while that for the total pressure falls in the range of 18-24 atm. The algebraic constraints contain X h , which denotes the vector of all process states, and U h , which denotes the vector of all process inputs. The

(23)

For a nonlinear MPC, the predicted output y over the prediction horizon P is obtained by the numerical integration of

(20)

subject to

(22)

subject to

φ ) ||PM1 - P*M1||2 + P*t ||2

∑||Γ[y(tk+i) - R(tk+i)]||2 +

∆u(tk),...,∆u(tk+M-1)i)1

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

(24)

y ) g(x)

(25)

from tk up to tk+P, where x and y represent the states and the output of the model, respectively. The symbol ||‚|| denotes the Euclidean norm, k is the sampling instant, Γ and Λ are diagonal weight matrices, and R ) [r(tk+1), ..., r(tk+P)]T is a vector of the desired output trajectory. ∆U(tk) ) [∆u(tk), ..., ∆u(tk+M-1)]T is a vector of M future changes of the MV vector ∆u that are to be determined by the online optimization. The control horizon (M) and the prediction horizon (P) are used to adjust the speed of the response and hence to stabilize the feedback behavior. Γ is usually used for tradeoff between different controlled outputs. The input move suppression, Λ, on the other hand, is used to penalize different inputs and thus to stabilize the feedback response. The objective function (eq 22) is solved online

Table 2. Process Parameters Bw Cp,p E Vg Tref

70 tonne 0.85 cal/s‚K 9000 cal/mol 500 m3 360 K

∆Hr kd MrCp,r Fg UA

-894 cal/g 0 s-1 1400 kcal/K 8500 mol/s 1.263 × 105 cal/s‚K

kp1 kp2 Tf ∆P ac

85 L/mol‚s 3 L/mol‚s 293 K 3 atm 0.548 mol/kg

Cp,H Cp,M1 Cp,N Cp,M2 Cp,w

7.7 cal/mol‚K 11 cal/mol‚K 6.9 cal/mol‚K 24 cal/mol‚K 1 cal/g‚K

2352


to determine the optimum value of ∆U(tk). Only the current value of ∆u, which is the first element of ∆U(tk), is implemented on the plant. At the next sampling instant, the whole procedure is repeated. To compensate for modeling error and eliminate steady-state offset, a regular feedback is incorporated into the output predictions, y(tk+1), through an additive disturbance term. Therefore, the output prediction is corrected by adding to it the disturbance estimates. The latter is set equal to the difference between plant and model outputs at present time k as follows:

d(tk) ) yp(tk) - y(tk)

(26)

The disturbance estimate, d, is assumed to be constant over the prediction horizon because of the lack of an explicit means of predicting the disturbance. However, for severe modeling errors, or open-loop unstable processes, the regular feedback is not enough to improve the NLMCP response. Hence, state or parameter estimation is necessary to enhance the NLMPC performance in the face of model-plant mismatch. In this work, an extended Kalman filter (KF) will be incorporated to correct the model state and, thus, to address the robustness issue. More details on the integration of an extended KF with the NLMCP algorithm are discussed elsewhere.14 Utilization of the NLMPC with KF requires adjustment of an additional parameter, σ. In a typical KF design, the diagonal elements of covariance matrices (Q and R) that represent the stochastic processes for unmeasured disturbance and measurement noise, respectively, are used to tune the filter. In the absence of accurate knowledge about the disturbance and noise characteristic, we further simplify the tuning method. Assuming Q ) q2I and R ) r2I and defining σ2 ) q2/r2 with r2 ) 1 will uniquely determine the KF gain and simplify its tuning to determining only one parameter.16 Intuitively, σ is the ratio between statistical measurement of the state uncertainty (disturbance) and the measurement noise. A large value for σ implies that the model uncertainty dominates over measurement noise. Therefore, σ is used to provide stability and robustness in the presence of modeling error, disturbances, and measurement noise. In addition to state estimation by KF, the predicted output will also be corrected by the additive disturbance estimates of eq 26. Online Linear MPC Algorithm Here the well-established linear model predictive controllers based on FIR models are used. A usual MPC formulation (notation follows that of Lee and Yu17) solves the following online optimization:

min

∆u(k),...,∆u(k+M+1)

||Γ[YP(k+1) - R(k+1)]||2 + ||Λ∆U(k)||2 (27)

subject to

YP(k+1) ) MPY(k) + SPm∆U(k) T

A ∆U(k) g b

the output vector, assuming ∆u(k+i) ) 0; i g M. Y(k) ) [y(k), ..., y(k+n-1)]T includes the predicted outputs over the truncation horizon n (length of FIR) based on ∆u(k+i) ) 0; i g 0. Equation 28 represents the output prediction based on the process model. The disturbance estimate defined by eq 26 should also be added in eq 28, or alternatively it can be absorbed in R(k+1). The latter is assumed for simplicity. As mentioned before, the disturbance is assumed to be constant over the prediction horizon. MP is a constant matrix consisting of ones and zeros, and SPm is what is usually referred to as the dynamic matrix of step response coefficients. The LMPC tuning parameters are the control horizon, M, the prediction horizon, P, the output weight, Γ, and the input weight (input-move suppression), Λ. Note that Λ can also be used to remove ill-conditioning. The regular feedback, i.e., additive disturbance, is valid only if the disturbances entering the system affect the output as independent random steps. For highly coupled multivariable systems or highly nonlinear systems with disturbances that tend to be slow drifts, LMPC with a constant (time-invariant) additive correction may fail to provide a reasonable performance. The performance is even worse for open-loop unstable processes. The well-established parameter/state estimation methods are based on state-space models. These techniques cannot be implemented directly in the employed LMPC formulation because the latter is based on FIR models. The alternative is thus to improve the model prediction by means of time-variant disturbance/modeling error predictions. For this reason, the LMPC robustness in this paper is improved by using a dynamic disturbance prediction model. The proposed method is based on input/output formulation.18 The following time series model for disturbance prediction

d(k) ) θ1d(k-1) + θ2d(k-2) + ... + θRd(k-R) (30) will be used in this paper. θ’s are the disturbance model coefficients to be identified online from actual modelplant mismatch estimates, and R is the size of the disturbance model, which is a parameter to be defined by the user. It acts as a forgetting factor and is used as a tuning parameter here. The larger R is, the more the disturbance model depends on the old values of the mismatch, d. In this paper, d(k) is, as defined before, the difference between the plant measurement and the model prediction at sampling instant k. The model coefficients can be determined directly by linear regression, i.e., solving a least-squares problem, from historical data for the mismatch, d. At each sampling instant and after θ’s have been updated, the disturbance model can be used to predict the future disturbance estimates via the following recursive equation:

(28) (29)

The symbols k, Γ, Λ, R(k+1), and ∆U(k) were defined before. YP(k+1) ) [y(k+1), ..., y(k+P)]T includes the predicted outputs over the future horizon P, where y is

R

d(k+l/k) )

θid(k+l-i) ∑ i)1

(31)

where the notation d(k+l/k) denotes the prediction at time k + l based on the actual measurement up to time

Ind. Eng. Chem. Res., Vol. 42, No. 11, 2003 2353

k, and the forecast estimates of the disturbances d(k+li) are defined as follows:

d(k+l-i) ) actual estimates of disturbance if l - i e 0 predicted disturbance based on measurement up to k if l - i > 0

{

}

The disturbance estimates are then added to the openloop output prediction, in addition to the additive disturbance, to minimize the model-plant mismatch and consequently enhance the LMPC performance. The new method requires adjustment of an additional parameter, R. A similar approach was presented by Sistu et al.19 and Ohshima et al.20 It should be noted that a disturbance prediction model is developed for each controlled output. Then, each disturbance model is used to correct the specific output for which the disturbance model was developed. It should be noted that for step change disturbance and at steady state, all of the past d(k) signals may take the same value. In this case, the determination of the model coefficient θ in eq 30 may suffer from singularity problems. When singularity occurs, the model coefficients are not updated. Instead, the last values for the model coefficients will be used. Polyethylene Process Control Loops Like any plant, the ethylene process contains several basic control loops. For example, the recycle flow is induced by the compressor and its flow rate is well controlled by a local controller that manipulates a control valve. Besides the low-level control loops, flow controllers, there are other important control loops that are directly related to product quality and process economics as discussed in the following. The total pressure and bed temperature controls are a common industrial practice because they preserve the total mass and energy balances and provide steady-state operation. Control of the partial pressures is essential because the partial pressure of the gases specifies the polymer grade. This multivariable control problem is a very challenging one because the gas partial pressures interact with each other. Also, there is a strong interaction between the partial pressure and the total pressure and between the monomer concentration and the bed temperature. Moreover, the temperature, total pressure, and flow rates are measured more frequently than the composition. The recycle stream has a large flow rate compared to other streams. Thus, it has a strong and fast influence on the process dynamics. The fresh feeds have, on the other hand, smaller flow rates and hence a slower effect on the dynamics, while the bleed flow has a fast dynamic on the total pressure but a slow dynamic on the partial pressure. Overall Control Schemes. It is known that a recycle stream creates significant control difficulty in the process. In fact, the recycle stream induces feedback interaction and increases the overall time constant of the process. It causes time-scale multiplicity and may even create instability. This behavior is known usually to occur for a recycle that comes from another processing unit in the plant and affects the mass balance directly. Although this is not exactly the same situation in this process, the differential time-scale phenomenon is also observed here. The different time scale occurs in this process because of the reaction kinetics, gas-phase

concentration, and feed rates. In fact, the time constant for the process variables ranges from 1.2 to 15 h. Kumar and Daoutidis12 addressed the control problem of a process with two time scales by decomposing the process dynamics into separate fast and slow time scales. The speed of process dynamics may segregate the process control loops into two groups; fast and slow. The partial pressure control loops belong to the second group, while the bed temperature and total pressure belong to the first group. In this work, we compare two overall control schemes. (1) A fully multivariable system that covers all of the six controlled variables is considered to be the first scheme and denoted as S1. In this case, the MVs are the feed flow rates, the bleed flow, and the coolant inlet temperature. (2) Another scheme, defined as S2, considers the partial pressures in a multivariable system and both the bed temperature and total pressure in separate external SISO loops. The MVs are the feed flow rates for the different partial pressures in a multivariable framework, the bleed flow for the total pressure in a SISO loop, and the coolant inlet temperature for the bed temperature in a SISO loop. The multivariable control system is conducted using NLMPC with a sampling rate of 20 min, while the individual SISO loops are handled by a PI controller with a sampling rate of 4 min. It should be noted that the catalyst feed rate has a strong effect on the reaction activity and hence on the mass and enthalpy balances of the reactor. Therefore, it can be used as MV. However, this is avoided in this work because, from a practical point of view, the catalyst feed system is very difficult to manipulate or, in the best case, imprecisely manipulated. It should be noted that, in practice, gas composition measurements may have time delays. However, the effect of time delays is not included in our analysis. Closed-Loop Simulations The polyethylene reactor operation is liable to several sources of disturbances. For simulation purposes, we classify the origin of the disturbances into three sources. One source is related to the catalyst feeding system and activity inside the reactor. Another source is the pressure system inside the reactor, which can be related to the recycle gas pump around the system, to periodic polyethylene removal, and to the bleed flow system. A third source is the thermal system, which can be affected by the ambient temperature, the temperature of the fresh gas feeds, and the external cooling system. For the first category, the effect of disturbances in the catalyst feed flow will be investigated. For the third category, the influence of upsets in the gas feed temperature and in the cooling water flow rate will be tested. No disturbance for the second category will be examined in this paper. Most of these disturbances are measured. However, for demonstration purposes they are assumed to present unmeasured sources of disturbances that affect the reactor operation. Consequently, their effect on the process dynamics and the effectiveness of the control system to improve the reactor operation will be examined. The following simulations, i.e., Figures 2-4, demonstrate the effect of the aforementioned disturbances on the process operation using scheme S2. In these figures, the dash-dotted lines show the process response when the partial pressures are not controlled. The dashed

2354


Figure 2. Feedback response to a -10% change in the coolant flow rate using scheme S2. Dash-dotted line: without partial pressure control. Dashed line: with partial pressure control.

Figure 3. Feedback response to a +3 °C change in the temperature of the fresh feed using scheme S2. Dash-dotted line: without partial pressure control. Dashed line: with partial pressure control.

lines show the response when the NLMPC is used to control the partial pressures. Note that, in both cases, the total pressure and the bed temperature are always controlled via independent PI control loops. The purpose of the simulations is to illustrate the effect of the disturbances on the process and to assess whether a controller for the partial pressure is needed or not. The parameters of the total pressure control loop are determined by trial and error to be kc ) 10 (mol/s)/atm and τI ) 0.1 s, while those for the temperature control loop are kc ) 1 °C/°C and τI ) 0.05 s. The NLMPC param-

eters are M ) 1, P ) 2, Λ ) [0] and Γ ) [1, 1, 1, 0, 1]. Note that the weight on the nitrogen partial pressure is taken as 0, which leaves PN without regulation. Despite this fact, the nitrogen partial pressure will eventually return to its target value. This is because the total pressure and the other partial pressures are controlled. In fact, the set point for PN has no significance in the real practice. The following arbitrary constraints on the MVs are used: 0 e FM1 e 260, 0 e FM2 e 7, 0 e FH e 3, 0 e FN e 5, 5 e Twi e 40, and 0.1 e Bt e 30. It should be noted that all of the following NLMPC simulations are based on small values for P (i.e., 2). This is because it is found through side simulations that larger values for P do not provide significant improvements. Figure 2 shows the process performance for a step change in the cooling water flow rate of -54 000 g/s, which corresponds to a 10% decrease in the flow. This upset is expected to reduce the cooling rate, leading to an increase in the reactor temperature. The resulting temperature variation will, in turn, affect the reaction activity and, consequently, the mass balance of the gas species. However, this amount of upset in the cooling rate led to a minor effect, as shown by the dash-dotted lines. In fact, the open-loop response of the partial pressures suffered from a minor overshoot. The dashed line shows how the NLMPC improved the reactor operation through reducing the overshoot in the controlled partial pressures, i.e., PM1, PM2, and PH. Of course, closing the partial pressure control loop brings its effect on T, Pt, and the uncontrolled PN. Figure 3 shows the process performance for a step change in the temperature of the gas fresh feeds of +3 °C, which is a common upset as the feed temperature fluctuates with ambient temperature. This type of


Figure 4. Feedback response to a +10% change in the catalyst flow rate using scheme S2. Dash-dotted line: without partial pressure control. Dashed line: with partial pressure control.

disturbance is expected to affect the enthalpy balance of the reactor. It turned out that this type of upset has a trivial influence on the process operation even in openloop mode, as shown by the dashed and dash-dotted lines. The reason for the marginal effect is that the thermal capacitance of the fresh feed is small compared to that of the recycle because of the large ratio of their corresponding flow rates. The dashed lines show how PM1, PM2, and PH are perfectly regulated when the NLMPC is implemented. The large scale in Figure 3 is used for a fair comparison with the other disturbance cases. Figure 4 illustrates the impact of an upset in the catalyst flow rate of 0.2 kg/s on the process dynamics. This amount of upset corresponds to a 10% increase in the catalyst flow rate, and it is expected to alter the reaction activity. Clearly, the increase of the catalyst flow rate had a considerable effect on the partial pressure as shown by the dash-dotted lines. The influence of the disturbance on the bed temperature and the total pressure is minor because they are controlled via the PI loops. In fact, a higher catalyst feed rate leads to a higher consumption of the ethylene, diminishing its partial pressure. Conservation of the total pressure causes the partial pressure of the other gases to build up. This situation will divert the polymer grade from its required value. The NLMPC is applied to overcome this problem. The dashed lines show that the detrimental effect of the disturbance is eliminated by the use of the NLMPC controller. Therefore, a smooth and consistent polymer grade can be maintained. In all of the above simulations, the total pressure and bed temperature were always maintained at their set point over the given range of operating conditions using the PI controller. According to the simulation results of Figure 4, the disturbance in the catalyst feed rate leads to an increase in the polymer production rate but

at the expense of a higher ethylene feed rate. In fact, the yield of polymer per ethylene fed is almost the same before and after the partial pressure control. This situation does not necessarily provide economical benefits. Comparison of the regulatory responses to the three different upsets highlights the following points. Upsets in Fw and Tf created a smaller effect on the process than upsets in Fc did. It can be seen that disturbances in Tf delivered very a marginal effect, and the reason for that was discussed earlier. More important is the open-loop response of the partial pressures for upsets in Fw and Tf. Interestingly, the partial pressure of the gases returns to their target values even though the NLMPC is not in service. These upsets impact the bed temperature directly and the partial pressure indirectly through the former. Because the bed temperature is well controlled and so is the total pressure, the partial pressures will automatically return to their steady-state values. On the other hand, the upset in the catalyst feed rate made a significant influence on the partial pressure even though the bed temperature and total pressure were perfectly controlled. This behavior is attributed to the fact that Fc has a direct affect on the partial pressure through the reaction activity. Therefore, in the case of disturbances that affect the partial pressure directly and not through the bed temperature, implementation of a control system for the partial pressure is a must to maintain a consistent product quality. This is true as long as the process is operating at a stable operating condition for the bed temperature. It remains to investigate which control scheme, i.e., S1 or S2, might provide the best control performance for this process. The following simulations compare the two schemes for the same regulatory control problems discussed above. The two schemes are compared using the same tuning parameter values of their controllers

2356


Figure 5. Feedback response to a -10% change in the coolant flow rate. Dashed line: scheme S1. Solid line: scheme S2.

Figure 6. Feedback response to a -10% change in the catalyst flow rate. Dashed line: scheme S1. Solid line: scheme S2.

for a fair comparison. The effect of Tf is ignored simply for its triviality. Figures 5 and 6 illustrate the process closed-loop performance for step change disturbances in Fw of

-54 000 g/s and in Fc of +0.2 kg/s, respectively. The dashed line corresponds to scheme S1 and the solid line to S2. Both schemes were able to reject the influence of the disturbance and maintain the control objective.


Figure 7. Feedback response to a +10% change in the catalyst flow rate in the presence of moderate modeling error. Solid line: scheme S2. Dashed line: scheme S1.

However, it can be seen that the performance of S2 outweighs that of S1. In fact, the performance of S1 exhibits higher overshoot and longer oscillation especially in the case of disturbances in Fw, More importantly, the S1 response suffers in all cases from a prolonged settling time, which appears like an offset because a relatively short simulation time is used. The slow response is due to the unified sampling rate used for all control loops. Note that a sampling rate of 20 min is used in the fully multivariable scheme S1. Further tuning of the NLMPC parameters, which should be avoided to ensure a fair comparison, will improve the response speed marginally. Our investigation indicated that a faster response for scheme S1 can be obtained only at a faster sampling rate. However, the smallest admissible sampling rate is limited by the frequency at which the measurements of the partial pressures are available. The latter is sampled at 10 min in a typical plant. It seems that isolating the fast loops from the slow loops, as suggested by Kumar and Daoutidis,12 gives a better solution to this two-time-scale control problem provided that the fast loops use a small sampling time. In fact, a much faster sampling rate is also possible for these loops because the measurement of the bed temperature and the total pressure is available at high frequency. It can be concluded that the performance of scheme S2 is better than that of S1. However, improvement of the S1 performance in terms of a faster response is possible if a smaller value for the sampling rate is allowable. Modeling Error. In addition to the above common disturbances, the effect of parametric uncertainty in the model will also be examined. For example, uncertainty in the reaction rate constants, kp1 and kp2, and in the activation energy, E, will be considered. These param-

eters are usually unknown exactly and have a strong effect on the polymerization reaction activity. Moreover, uncertainty in the catalyst active site concentration, ac, will also be investigated. The latter is also imperfectly known and influences the catalyst activity inside the reactor. Finally, uncertainty in the overall heat-transfer coefficient, U, of the external cooler will also be tested. This parameter has a strong impact on the cooler performance and hence can seriously affect the reactor temperature. Figure 7 illustrates the closed-loop response for a step change in the catalyst flow rate of 0.2 kg/s in the presence of modeling errors without incorporating KF. Specifically, the overall heat-transfer coefficient of the model, U, is considered to be 10% higher than that of the plant, the reaction rate constants, kp1 and kp2, of the model are taken as 10% lower than those of the plant, and the catalyst deactivation constant, kd, is taken to be 1 × 10-5 s-1 in the model instead of 0 (no deactivation) in the plant. The solid line represents the response of scheme S2, and the dashed line represents the response of scheme S1. When an imperfect model is used, the closed-loop performance became considerably inferior compared to that shown in Figure 4 by the dashed lines. The transient response of most outputs suffered from substantial overshoot, damped oscillation, and long settling time. Nevertheless, the NLMPC was still able to stabilize the process and to eventually maintain the controlled variables at their set points. A similar observation can be made for the closed-loop response to the other sources of disturbances. Note that the controller was able to achieve its objective in the presence of modeling error even without using KF except for regular feedback corrections. The performance can be obviously improved substantially by KF. How-

2358


Figure 8. Feedback response to a +10% change in the catalyst flow rate in the presence of a moderate modeling error plus a +50% change in E using scheme S2. Solid line: with KF. Dashed line: without KF.

ever, we choose not to use the filter to demonstrate the ability of the controller to handle this degree of modeling error in the model. As the dashed line demonstrates, scheme S1 is still inferior to S2 in the sense of slower settling time. Next we investigate the effectiveness of incorporating KF into the NLMPC performance. To illustrate how KF can expand the effectiveness of NLMPC to handle a wider range of uncertainty, more severe uncertainty conditions will be examined. The following simulations examine the NLMPC closed-loop performance for an upset in Fc of +0.2 kg/s in the presence of the same moderate modeling errors as those used in Figure 7 plus an additional one or a larger magnitude of error in a specific process parameter, as described in each of the following simulations. Only scheme S2 will be examined. The purpose of the KF is to update all of the model states using the plant measurement of the gas compositions and reactor temperature. The KF is tuned by adjusting its parameter σ to 1. This value is found by trial and error and is fixed in all simulations. Experience indicates that σ should be of order unity or less and that when it is around the optimal value, σ can vary by at least an order of magnitude with little change in performance.21 Figure 8 examines the above regulatory problem augmented by additional uncertainty in the activation energy, E. The activation energy of the model is taken to be 50% higher than the nominal value used in the plant. This amount of parametric error created a clear model-plant mismatch, causing an offset and persistence oscillations when filtering is not employed, as demonstrated by the dashed line. When the model states are corrected properly, NLMCP succeeded in stabilizing the process with some initial loss in the performance,

as shown by the solid lines. It is clear that the bed temperature and total pressure can be easily brought to steady-state operation. However, loss of the polymer grade might be observed for at least 10 h because of the oscillations in the partial pressures. The effect of -50% inaccuracy in the model activation energy is also tested. Good control was obtained using KF; however, the simulation is excluded because the feedback response without filtering was not as bad as that of the previous case. Figure 9 tests the MPC performance for the moderate modeling errors with an elevated uncertainty on kd. Specifically, the catalyst deactivation in the model made 1 × 10-3 s-1. Increasing the catalyst feed rate to the plant and increasing the catalyst deactivation in the model widen the discrepancy in the number of active catalyst moles between the plant and model. If no KF is used, the NLMCP feedback response becomes unstable in terms of obvious excursions, as shown by the dashed line in Figure 9. For the same NLMPC parameter values and with the aid of KF, the feedback response is stabilized and the performance is enhanced, as shown by the solid line in the same figure. The temperature and total pressure returned rapidly to their target values because they have fast dynamics. The partial pressures of PH, PN, and PM2 required a longer time to settle back to their target values. In the same way, the impact of (50% uncertainty in the values of both kp1 and kp2 on the control performance is examined. The error in the other parameters is the same as that in Figure 7. The simulation indicated that the detrimental effect of the modeling errors intensifies at -50%. A similar conclusion about the strength of model corrections by KF can be drawn. The results are not shown because no significant behavior over the


Figure 9. Feedback response to a +10% change in the catalyst flow rate in the presence of a moderate modeling error but at kd ) 10-3 using scheme S2. Solid line: with KF. Dashed line: without KF.

previous results was detected. Similarly, the effect of an additional error of (50% on the catalyst active concentration, ac, was examined. However, insignificant results were obtained. Therefore, the simulations are omitted. The controllers is also examined for a (50% error in the model overall heat-transfer coefficient, U, while the error in the other parameters is the same as that in Figure 7. The simulation in Figure 10 is for the case of +50% error in U. The dashed line in Figure 10 depicts the NLMPC response when the controller relies only on the regular feedback correction. As a result, a loss of control is observed and is manifested in the continuous oscillation. On the other hand, the solid line in Figure 10 attests to the effectiveness of KF in improving the NLMPC performance. Some minor overshoot occurs initially, but the controller was able to bring the process to steady state within 10 h. For the case of -50% deficiency in U, the process becomes very difficult to control, as shown in Figure 11. The dashed line in Figure 11 corresponds to NLMPC without KF. The control performance is somewhat poor especially for PH and PN. Although the responses might seem acceptable, they deteriorate for longer simulation time. Note also that the simulation works only at P ) 1. Other values for P will make the simulation crash because of instability. The solid lines in Figure 11 correspond to NLMPC with KF. Obviously, there is no recorded improvement. The simulation is carried out at P ) 2. At other values for P, the simulation either is inferior or crashes as a result of instability. Our investigation revealed that the controller may derive the bleed flow to zero or to a very small value to counter the effect of the disturbance. If the bleed flow is kept at a small value for a long period, the process behaves like an integrator system, leading

to a quick buildup of the total pressure. As a result, the simulation may crash. The last test presents a challenging control problem. At the plant, the elevated catalyst flow rate increases the bed temperature and reaction activity substantially, leading to a quick consumption of the ethylene. In the model, the situation is even more serious because the cooler capacity is limited as a result of the reduction of the overall heat-transfer coefficient. This limitation leads to escalation in the bed temperature, which will make the temperature control loop bring Twi to its lower limit. The high bed temperature and high catalyst feed rate will lead eventually to excessive consumption of the ethylene. As a result, the total pressure falls and the bleed flow goes to zero to compensate for the pressure drop. The dark solid line in Figure 11 shows the process response for the same test but at a sampling time of 10 min for the NLMCP algorithm with KF. The sampling time is reduced for the sake of a better controller performance. Although the performance is improved substantially, the partial pressure of M1, H2, and N2 was not brought completely to its set point. In fact, it is found that further speeding of the NLMPC execution, exactly at the same rate at which the PI is implemented, introduces disturbances to the PI control loops. Our investigations revealed that further control improvement is not possible because the PI control loops are isolated from the NLMPC loop. This means that the entire available MVs are not utilized optimally or collectively. For example, when Twi and/or Bt of the PI loops are stuck at their lower limits, the other MVs were not sought for possible improvement of the PI control loops. This can be achieved only if a fully multivariable control system is employed. Figure 12 demonstrates the

2360


Figure 10. Feedback response to a +10% change in the catalyst flow rate in the presence of a moderate modeling error but at a +50% change in U using scheme S2. Solid line: with KF. Dashed line: without KF.

Figure 11. Feedback response to a +10% change in the catalyst flow rate in the presence of a moderate modeling error but at a -50% change in U using scheme S2. Solid line: with KF. Dashed line: without KF. Dark solid line: with KF and a sampling time of 10 min.

performance of the fully multivariable control scheme (S1) with a unified sampling time of 10 min. It is obvious that the controller successfully improved the overall

performance. Moreover, the performance became superior to that of S2 shown in Figure 11. The initial response suffers from some oscillation and drift from


Figure 12. Feedback response to a +10% change in the catalyst flow rate in the presence of a moderate modeling error but at a -50% change in U. Solid line: using scheme S1 with KF and Ts of 10 min. Dashed line: using scheme S1 without KF and Ts of 20 min.

set point, which may create an off-specs product for about 10-20 h. However, this good control was achieved at extreme values for the modeling errors. Optimal Operation. The discussion of the tests shown in Figures 2-4 indicated, as mentioned earlier, the need for optimizing the bleed flow. The discussion of the results shown in Figures 5 and 6 pointed out that the performance of scheme S1 can be made as good as that of S2 if a smaller sampling time is used. Finally, the results of Figure 12 imply the superiority of S1 over S2 in the case of disturbances that drive the MVs to the edge of constraints. Therefore, it might be interesting to utilize the full multivariable NLMPC (scheme S1) to optimize the ethylene conversion by operating the plant at the same grade but at a lower bleed flow, thus presenting one more advantage of S1. One advantage of the multivariable NLMPC algorithm is the ability to operate the process around optimal conditions. For example, the loss of the raw ethylene or equivalently the ethylene conversion can be improved. One way to improve the ethylene conversion is to cut down the bleed flow. Therefore, we test scheme S1 for maintaining the process variables at their desired targets while reducing the bleed flow to another smaller steady-state value. In this case, the set point for the bleed flow will be included as a controlled variable and its valve position as the MV. The result of the simulation is shown in Figure 13. The solid line in the figure shows the process response to a -50% set point change that occurs after 140 min from the beginning of the simulation. The simulation also includes the moderate modeling errors of Figure 7. Because the valve dynamics of the bleed flow is not modeled here, the response of the bleed flow is instantaneous as shown in the figure. This sudden change in the bleed flow creates disturbances

Figure 13. Set-point change of -50% in the bleed flow using scheme S1 in the presence of a moderate modeling error.

to the other variables. As a result, the controller regulates the other MVs simultaneously to maintain the control objectives. Therefore, some transient is observed in controlled variable responses. A slow drift is also detected in the PH and PN responses; however, the drift is reasonable in magnitude and disappears after some time. The regulation of FM1 is shown in the figure, while that for the other MVs is omitted for brevity. To assess the benefit, the calculated conversion at the initial steady state is 96.56% and that for 50% reduction in

2362


Figure 14. Feedback response to a +10% change in the catalyst flow rate using LMPC with scheme S2.

the bleed flow is 98.18%. This indicates that a 1.6% increase in the conversion can be obtained. Linear MPC Tests. So far, only NLMPC was examined. However, LMPC is more preferable industrially. Therefore, in this section we test the LMPC for the same control problem. Specifically, the effectiveness of LMPC for the rejection of 0.2 kg/s upset in the catalyst flow rate using schemes S1 and S2 will be investigated. This type of disturbance is chosen because it presented the most difficult case. There are some issues related to the successful implementation of LMPC. First, the size of the FIR model is crucial. Because the process has slow dynamics, then a large truncation number is needed to represent the steady-state behavior when a small sampling time is used. Therefore, the truncation number imposes an additional limitation on the smallest possible sampling rate. Moreover, the linearized process model is ill-conditioned, which requires more attention to be paid for the tuning of the LMPC parameters, as will be shown in the simulations. Although the FIR model is derived from a perfect process model, i.e., no parametric modeling errors are considered, some degrees of model-plant mismatch exist because of linearization. Therefore, special attention should be made for improving the model predictions, as will be demonstrated. Figure 14 illustrates the LMPC performance when scheme S2 is implemented. For this case, an FIR model is developed from the nominal plant using a sampling time of 20 min and a truncation number of 120. Because we are testing scheme S2, the model was developed from the nominal plant under T-Twi and Pt-Bt control loops. The figure shows stable feedback responses. However, the response of PM1, PH, and PN suffers from apparent steady-state offset. The fast dynamics loop, i.e., T and Pt, which relies on the isolated PI algorithm delivered

offset-free responses. It should be noted that the overall closed-loop response was stabilized only by both careful tuning and incorporation of the proposed disturbance estimation procedure. In fact, it was necessary to adjust the move suppression matrix Λ to diag[10, 20, 20, 30] in order to eliminate the ill-conditioning. In addition, a large prediction horizon, P, was also essential to provide a stable and fast response. Specifically, P ) 100 was used. Values for P of less than 50 were found to provide an unacceptable sluggish performance, while values larger than 100 delivered a very aggressive response. Moreover, it was found that relying on the regular feedback correction delivered closed-loop responses with a considerable steady-state offset. The offset is due to model-plant mismatch because the plant is nonlinear and the model is a linear FIR model. Incorporation of the proposed disturbance estimation technique with R ) 5 helped in substantially reducing the offset. The above values for the tuning parameters and R were found by an intensive stepwise trail-and-error procedure. Our investigation revealed that further improvement of the closed-loop performance is not possible. Any attempt to adjust the LMPC parameter values and/or R will result in unstable responses. Note that the previous NLMPC simulations did not require tuning of its parameters for stability. In fact, NLMPC was found to deliver stable responses for any value for the tuning parameters even for the most aggressive ones, i.e., P ) M ) 1 and Λ ) [0]. Figure 15 illustrates the scheme S1 performance for the same control objective. For this case, an FIR model was developed from the nominal plant using a sampling time of 10 min and a truncation number of 240. The FIR model is developed from a fully open system; i.e., no control loops are in service. The sampling time was intentionally made shorter to improve the performance.


Figure 15. Feedback response to a +10% change in the catalyst flow rate using LMPC with scheme S1.

Even so, the control problem became more challenging than that in Figure 14. Specifically, the linear model that covers the six controlled variables is highly illconditioned. In fact, the condition number is increased by manyfolds from that of the model used in S2. It was extremely difficult to find a set of values for Λ that eliminate the ill-conditioning by trial and error. Alternatively, the optimum values for the diagonal elements of Λ that reduce the ill-conditioning are found by solving an objective function that minimizes the condition number as follows:

min |cond(Kmpc) - 1| ) |cond[(ΓTΓ +

diag(Λ)

SPmTΛTΛSPm)-1SPmTΓ] - 1| In the above objective function, Kmpc is the unconstrained LMPC controller gain.22 The above objective function is solved using available optimization software to find the diagonal elements of Λ that brings the condition number of Kmpc as close as possible to one at Γ ) I. The above problem is solved at different values for M and P, which indicated that the condition number decreases with increasing P and decreasing M. The optimal value was found to be Λ ) [0.7178, 0.9740, 0.6546, 0.5516, 4.6747, 0.550] at M ) 1 and P ) 100. These tuning parameter values were enough to provide stable responses but not good performances. The next step was to adjust Γ so that a better performance, i.e., less offset, can be obtained. Interestingly, a strong tradeoff between the controlled variables prevails, as shown by the dashed and solid lines in Figure 15. For example, when Γ is tuned to reduce the offset in PM1, the offset in the other variables, especially T and Pt, is increased and vice versa. The best achievable performance is shown by the dash-dotted lines in Figure 15, which was obtained at Γ ) [4.5, 1, 10, 0, 20, 10]. Note that the disturbance estimation model was also incorporated with slight modification. In the previous simulation (Figure 14), the disturbance estimation was implemented with a unified R for all outputs. This was found to provide a marginal enhancement in the current simulation. Therefore, different values for R were used for each disturbance model. This modification managed to help improve the closed-loop performance. Specifically, R for each disturbance model was tuned to [3, 30,

5, 0, 1, 0]. The necessity for variable values for R is attributed to the uneven response time of the controlled outputs. Despite all attempts, the resulting feedback response shown in Figure 15 still suffers from steadystate offset. Note that the simulation time in Figure 15 was made deliberately long to clearly illustrate the existence of static offset. It was shown in the last two simulations that incorporating the proposed disturbance prediction is a must to improve the performance and in some cases to stabilize the response. This situation is observed even though the FIR model was developed for the nominal plant. If the model was developed from an imperfect plant model or the plant dynamics deviates form the nominal plant dynamics, the LMPC performance may even worsen. It is also clear that implementing LMPC in scheme S1 is inferior to scheme S2, which coincides with the finding that was obtained from the NLMPC simulations. The inferiority of S1 compared to S2 is evident from the larger offset shown in Figure 15 and from the difficulty in tuning and implementing scheme S1. Conclusions This paper addresses the multivariable control problem of an industrial polyethylene reactor. The control problem is difficult because of the differential response times of the process controlled variables. In addition, these variables are measured at a differential sampling rate depending on the nature of their sensors. The latter places a restriction on the smallest possible sampling rate that can be used for control executions. For this reason two control schemes were examined and compared. One scheme handles the whole controlled variables by one multivariable control algorithm, LMPC or NLMPC. The other scheme segregates the controlled variables into a fast dynamic loop handled by a multiloop PI algorithm and a slow dynamic loop handled by a multivariable control algorithm, i.e., LMPC or NLMPC. The comparison favored the second scheme in terms of faster settling time and less overshoot in some cases. This was found to be true whether LMPC or NLMPC is used. However, there was a clear advantage for using the first scheme in a number of specific cases. One case is during the saturation of the MVs of the fast dynamic loop. Another case is during the transition from one operating condition to another optimal one. The comparison of LMPC with NLMPC indicated the superiority of the latter. NLMPC provided a better performance in the sense of offset-free responses. NLMPC implementation was easier because LMPC necessitated careful and intensive tuning. NLMPC was more robust in the face of plant-model mismatch. NLMPC provided a good performance using a simple regulatory feedback even in the presence of a moderate modeling error. When KF is incorporated, NLMPC was able to provide a good performance over a wide range of modeling errors. On the other hand, the performance of LMPC can only be improved, to some extent, by incorporating the disturbance estimation method. This was necessary even though the FIR model was developed from the perfect process model. In addition, the time-scale multiplicity of the process leads to an ill-conditioning problem, which was more paramount in the LMPC approach. For this reason, careful tuning of the move suppression matrix was a must to stabilize the LMPC performance. Therefore, we can say that a state estimation procedure

2364


should be used to improve the NLMCP performance while a disturbance/modeling error estimation must be used to improve the LMPC performance.

σ ) tuning parameter for the Kalman filter R ) tuning parameter for the disturbance prediction model θ ) coefficient of the disturbance prediction model

Nomenclature

Literature Cited

A ) constant matrix for linear constraints ac ) active site concentration, mol/kg b ) Vector of upper and lower bounds for the linear constraints Bw ) mass of the polymer in the bed, tonne Bt ) bleed flow rate, mol/s CM1, CM2, CN, CH ) concentrations of monomer, comonomer, nitrogen, and hydrogen, mol/m3 Cp,M1, Cp,M2, Cp,H, Cp,N ) heat capacities of monomer, comonomer, hydrogen, and nitrogen, cal/mol‚K Cp,g, Cp,w ) heat capacities of recycle gas and water, cal/ mol‚K Cp,p ) heat capacity of polymer, cal/g‚K E ) activation energy for propagation, cal/mol Fc ) catalyst flow rate, kg/s Fw, Fg ) cooling water and recycle flow rates, mol/s FM1, FM2, FN, FH ) monomer, comonomer, hydrogen, and nitrogen flow rates, mol/s HF, HG, HP ) enthalpy of fresh feed, recycle gas, and product, cal/s HT ) enthalpy of the gas leaving the reactor, cal/s HR ) enthalpy change due to polymerization, cal/s kd ) deactivation rate constant, s-1 kp1, kp2 ) propagation rate constant for monomer and comonomer, L/mol‚s M ) control horizon MP ) constant matrix Mw ) water holdup in the heat exchanger, mol MrCp,r ) thermal capacitance of the reaction vessel, kcal/K n ) truncation number for FIR models Op ) polymer outlet rate, kg/s P ) prediction horizon Pt ) total pressure, atm PM1, PM2, PN, PH ) partial pressures of monomer, comonomer, nitrogen, and hydrogen, atm R ) ideal gas constant, atm‚m3/K‚mol R ) vector of set points RM1, RM2, RH ) consumption rates of monomer, comonomer, and hydrogen, m3/mol‚s P Sm ) dynamic step response matrix T, Tf, Tref ) bed, feed, and reference temperatures, °C Tgi, Tg ) temperature of the recycle stream before and after cooling, °C Twi, Two ) cooling water temperature before and after cooling, °C t ) time, s UA ) overall heat-transfer coefficient multiplied by the heat-transfer area, cal/s‚K Vg ) gas holdup in the reactor, m3 x ) vector of states xM1, xM2, xN, xH ) mole fractions of monomer, comonomer, nitrogen, and hydrogen Y, YP ) vector of future outputs over n and P, respectively Yc ) number of moles of catalyst sites, mol y, yp ) vector of model outputs and of plant outputs

(1) Choi, K. Y.; Ray, W. H. The Dynamic Behavior of Fluidized Bed Reactors for Solid Catalyzed Gas-Phase Olefin Polymerization. Chem. Eng. Sci. 1985, 40, 2261-2279. (2) McAuley, K. B.; McDonald, D. A.; Mclellan, P. J. Effects of Operating Conditions on Stability of Gas-phase Polyethylene Reactors. AIChE J. 1995, 41, 868-879. (3) Gasem, N. Effect of Polymer Particle Size and Inlet Gas Temperature on Industrial Fluidized Bed Polyethylene Reactors. Chem. Eng. Technol. 2000, 22, 777-787. (4) Ajbar, A.; Alhumaizi, K. Gas-Phase Polyethylene Reactors a Global Study of Stability Behaviour. Chem. Eng. Res. Des. 2001, 79, part A, 195-208. (5) Dadebo, S.; Bell, M.; Mclellan, P.; McAuley, K. Temperature Control of Industrial Gas-Phase Polyethylene Reactors. J. Process Control 1997, 7, 83-93. (6) Ali, E.; Abasaeed, A.; Al-zahrani, S. Optimization and Control of Industrial Gas-Phase Ethylene Polymerization Reactors. Ind. Eng. Chem. Res. 1998, 37, 3414-3423. (7) Ali, E. M.; Abasaeed, A. E.; Al-Zahrani, S. M. Improved Regulatory Control of Industrial Gas-Phase Ethylene Polymerization Reactors. Ind. Eng. Chem. Res. 1999, 38 (6), 2383-2390. (8) Seki, H.; Ogawa, M.; Ohshima, M. PID Temperature Control of an Unstable Gas-Phase Polyolefin Reactor. J. Chem. Eng. Jpn. 2001, 34, 1415-1422. (9) McAuley, K. B.; MacGregor, J. Nonlinear Product Quality Control in Industrial Gas-phase Polyethylene Reactor. AIChE J. 1993, 39, 855-866. (10) Luyben, W.; Tyreus, B.; Luyben, M. Plant Wide Process Control; McGraw-Hill: New York, 1999. (11) Morud, J.; Skogestad, S. Effects of Recycle on Dynamics and Control of Chemical Processing Plants. Comput. Chem. Eng. 1994, 18, 529-534. (12) Kumar, A.; Daoutidis, P. Nonlinear Dynamics and Control of Process Systems with Recycle. J. Process Control 2003, in press. (13) Christofides, P. Output Feedback Control of Nonlinear Two-Time-Scale Processes. Ind. Eng. Chem. Res. 1998, 37, 18931909. (14) Ali, E.; Zafiriou, E. Optimization-based Tuning of Nonlinear Model Predictive Control with State Estimation. J. Process Control 1993, 3, 97-103. (15) McAuley, K. B.; MacGregor, J. F.; Hamielec, A. E. A Kinetic Model for Industrial Gas-Phase Ethylene Copolymerization. AIChE J. 1990, 36, 837-850. (16) Hamilton, J. C.; Seborg, D. E.; Fisher, D. G. An Experimental Evaluation of Kalman Filtering. AIChE J. 1973, 19, 901908. (17) Lee, J. H.; Yu, Z. H. Tuning of Model Predictive Controllers for Robust Performance. Comput. Chem. Eng. 1994, 18, 15-37. (18) Ali, E. Robust Model-Based Control of Open-loop Unstable Processes. Arabian J. Sci. Eng. 1999, 24 (1B), 33-48. (19) Sistu, P. B.; Gopinath, R. S.; Bequette, B. W. A Disturbance Estimator for Model Predictive Control. In Advanced Control of Chemical Processes; Bonvin, D., Ed.; Pergamon Press: New York, 1994; pp 173-178. (20) Ohshima, M.; Ohno, H.; Hashimoto, I.; Sasajima, M.; Maejima, M.; Tsuto, K.; Ogawa, T. Model Predictive Control With Adaptive Disturbance Prediction and its Application to Fatty Acid Distillation Column Control. J. Process Control 1995, 5, 41-48. (21) Ricker, N. L. Model Predictive Control with State Estimation. Ind. Eng. Chem. Res. 1990, 29, 374-382. (22) Lee, J. H.; Morari, M.; Garcia, C. E. State-Space Interpretation of Model Predictive Control. Automatica 1994, 30, 707717.

Greek Letters ∆u ) vector of manipulated variables ∆U ) vector of M-future manipulated variables ∆Hr ) heat of reaction, cal/g Λ ) input weight Γ ) output weight

Received for review May 23, 2002 Revised manuscript received February 10, 2003 Accepted March 19, 2003 IE020386H

Multivariable Control of a Simulated Industrial Gas-Phase

Recommend Documents