Ind. Eng. Chem. Res. 1997, 36, 389-398
389
Nonlinear Model Predictive Control of Industrial Type IV Fluid Catalytic Cracking (FCC) Units for Maximum Gasoline Yield Emad E. Ali* Department of Chemical Engineering, King Saud University, P.O. Box 800, Riyadh 11421, Saudi Arabia
Said S. E. H. Elnashaie† Chemical and Petroleum Engineering Department, United Arab Emirates University, P.O. Box 17555, Al-Ain, United Arab Emirates
This paper addresses the problem of stabilizing the operation of an industrial type IV fluid catalytic cracking unit (FCCU) around an unstable high gasoline yield operating steady state. The model used in the work has been previously checked successfully against a number of industrial units. The stabilization is achieved through simulated implementation of a model predictive control (MPC) version that utilizes the nonlinear model in the output prediction. A state estimation technique is also incorporated in the MPC algorithm to improve the latter performance in the presence of unmeasured disturbances and/or parametric modeling errors. Two types of closed-loop simulations were tested, namely, servo and regulator problems. The results of the simulation illustrated that operating the FCCU at the unstable region is possible due to the effectiveness of MPC in handling such a problem. It is also found that, for an openloop unstable process, state estimation must be used to improve the feedback response of MPC in the presence of unmeasured disturbances or modeling errors. Traditional PI control was also examined for the FCCU stabilization control problem, and its performance is compared with that of the NLMPC. The investigation revealed that reasonable PI performance can also be obtained through careful tuning of its settings and that NLMPC is considered easier to apply in that case. 1. Introduction Fluid catalytic cracking (FCC) is one of the most important units in the petroleum refining industry for the conversion of heavy gas oil to gasoline and light hydrocarbons. The performance of the FCC units plays a major role on the overall economics of petroleum refineries. Thus, any improvement in the operation or control of FCC units will result in considerable economic benefits. Several investigations have been concerned with the general behavior of the industrial FCC units. Some of these research efforts dealt with the modeling of FCC units which are useful in elucidating some of the main characteristics of these units for better design, operation, and control (Luyben and Lamb, 1963; Kurihara, 1967; Iscol, 1970; Lee and Kugelman, 1973; Elnashaie and El-Hennawi, 1979; Edwards and Kim, 1988; Elshishini and Elnashaie, 1990a,b). A type IV industrial FCC unit consists of two interconnected gas-solid fluidized beds. The first one is the reactor, in which gas oil is catalytically converted to gasoline, coke, and light hydrocarbon gases. The catalyst is regenerated in the second fluidized bed, the regenerator, by using air to burn the coke deposited on the catalyst in the reactor. The coke-burning step also supplies the heat for the endothermic cracking reactions through circulating the catalyst between the two vessels. The mathematical modeling of this unit is complicated by the high nonlinearity of the process. A rigorous mathematical model should consider the twophase nature of the fluidized beds in both the reactor * Corresponding author. E-mail:
[email protected]. Telephone: ++(9661)467-6871. Fax: ++(9661)467-8770. † On leave from the Chemical Engineering Department, Cairo University, Cairo, Egypt. E-mail:
[email protected]. Fax: ++(9713)645-083. S0888-5885(96)00357-0 CCC: $14.00
and the regenerator, together with their strong interactivity using reliable kinetics schemes. Different models for type IV FCC units of different degrees of sophistication were developed (Elnashaie and Elshishini, 1993a). Most of these mathematical models have shown that these industrial units are operating in the multiplicity region, specifically at the middle unstable steady state in order to achieve the maximum gasoline yield (Iscol, 1970; Elnashaie and El-Hannawi, 1979; Edwards and Kim, 1988). The static and dynamic bifurcation behaviors of the uncontrolled type IV FCC units have shown that this operating point is very sensitive to the variation of the process parameters (Elnashaie and Elshishini, 1993b). These findings suggest that a powerful control scheme should be used for the stabilization of this challenging control problem. Several earlier control studies on FCC units were reported (Kurihara, 1967; Prett and Gillette, 1980; Huang, 1989). Recently, there have been many efforts to investigate the FCC control problem through implementation of model predictive control (MPC) algorithms (Lokesh and Georgakis, 1994; Khandalekar and Riggs, 1993; Schmid and Biegler, 1990; Arbel et al., 1993; Cutler et al., 1992; Balchen et al., 1992; Grosdidier et al., 1993). Most of these efforts are based on suppressing the system response to external disturbances, while in fact the primary control problem associated with FCC units is the stabilization of the unstable operating state of the system. The aim of this paper is, thus, to present an investigation of the effectiveness of a nonlinear MPC (NLMPC) to operate these FCC units at their corresponding desirable unstable operating points. NLMPC has been used to stabilize the operation of unstable processes (Ali and Zafiriou, 1993; Hidalgo and Brosilow, 1990; Kent and Fisher, 1993). In the present work the investigation focuses on the implementation © 1997 American Chemical Society
390 Ind. Eng. Chem. Res., Vol. 36, No. 2, 1997
of the NLMPC version developed by Ali and Zafiriou (1993) to stabilize the unstable high gasoline yield operating steady states of industrial type IV FCC units. A state estimation technique is incorporated with NLMPC to improve its performance in the face of unmeasured disturbances and/or modeling errors. Other estimation techniques have been used by Wright et al. (1991), Li and Biegler (1990), Eaton and Rawlings (1990), and Jang et al. (1987) to estimate the unknown parameters of the model by using on-line optimization. However, in those investigations, the optimization index is solved for a specific time horizon which is determined by the user. This introduces a new tuning parameter and also increases the computational load if the prediction horizon is long. In addition, the degree of freedom of the optimization problem in these methods is the model unknown parameters which may pose a problem since it is difficult sometimes to determine or limit the number of unknown parameters. The proposed estimation technique, in the present work, differs from the one developed by Ali and Zafiriou (1993) in that computation of the steady-state filter gain is obtained via on-line optimization instead of iterative solution of the Riccati differential equation. This modification is used to avoid the possible failure in the computation of the filter gain due to loss of the observability condition or to the internal instability of the linearized model. The modified state estimation method keeps the same simple structure of the state estimation by kalman filtering and at the same time avoids the limitation associated with the model characteristics. As mentioned earlier, the modified estimation is based on on-line optimization, but it differs from the optimizationbased estimation techniques in a way that does not require the estimation of a large number of parameters nor intensive computations. Detailed description of the modified method is given in section 3. Simply, the proposed estimation method is a modified version of the one used earlier by one of the authors and combines the good features of two estimation techniques used in the literature. 2. Mathematical Model In this investigation a mathematical model is used which is based on the two-phase theory for bubbling fluidized beds. This model was developed earlier be Elnashaie and El-Hennawi (1979) and later modified and verified against industrial type IV FCC units by Elnashaie and Elshishini (1993b) and Elnashaie et al. (1995). The cracking reactions network used in this model is based on the three-pseudocomponent kinetic model of Weekman (1968). gas oil (A1)
K1
gasoline (A2)
K2
coke + gases (A3)
K3
The rates of reaction for this consecutive-parallel network are given by Weekman and Nace (1970). For the coke burning in the regenerator, the model uses the continuous reaction model for solids of unchanging size presented by Kunii and Levenspiel (1969). The FCC model used in this study recognizes the two-phase nature of the fluidized beds in the reactor and the regenerator with the assumption of pseudo steady state in the bubble phase, in addition to the other assumptions described earlier (Elnashaie and Elshishini, 1993b;
Table 1. FCCU Data fresh feed (kg/s) recycle HCO (kg/s) combined feed ratio (CFR) air flow rate (kg/s) combined feed temperature (K) air feed temperature (K) dimensionless air feed temperature (yfa) hydrogen in coke (wt %) HCO/(HCO + CSO + LCO) coke/(coke + gases) CO2/(CO2 + CO) (m3/m3) -∆Hc (kJ/kg) regenerator dimensions reactor dimensions catalyst retention in reactor (metric ton) catalytic retention in regenerator (metric ton) API of raw oil feed reactor pressure (kPa) regenerator pressure (kPa) average particle size (m) pore volume of catalyst (m3/kg) apparent bulk density (kg/m3) catalyst surface area (m2/kg)
16872 2.108 1.12565 10.67 527 436 0.872 4.17 0.286 0.221 0.75 31235.6 5.334 m i.d. × 14.859 m TT 3.048 m × 12.76 m TT 17.5 50 28.7 225.4938 224.8708 0.00072 0.00031 800 215
Elnashaie et al., 1995). The model includes the dynamics of gas oil and gasoline in the reactor, catalyst activity in the reactor, reactor dense phase temperature, catalyst activity, and dense phase temperature in the regenerator. The dynamic model equations are given in dimensionless form as follows:
dX1 ) B(Yv - X1) + a1(X2 - X1) + dt a2(Yf - Yv) - Hv + HRCB - HLR dX2 EHR ) [B (Y - X2) + a3(X1 - X2) + dt EHG D fa βc(CB) - HLG] dX3 1 ) [B(X3f - X3) - X5(R1e-γ1/X1 - R3e-γ3/X1)X32] dt Rs dX4 1 ) [B(X4f - X4) - X5(R2X4e-γ2/X1 - R1X12e-γ1/X1)] dt Rd dX5 EHR ) [C (X - X5) - X5Rc] dt EMR RD 6 dX6 EHR ) [C - CGD(X6 - X5)] dt EMG B The detailed derivation of the model is given elsewhere (Elnashaie and El-Hennawi, 1979; Elshishini and Elnashaie, 1990a,b). All variables and parameters are defined in these references and also in the Nomenclature section in the present paper. The static version of this model has been previously compared with two industrial catalytic cracking units, and it was shown that both units are operating at the middle unstable steady state and that the multiplicity region covers a very wide range of parameters (Elshishini and Elnashaie, 1990a,b). In the present paper, the design and operating data for one of these industrial units is used and given in Table 1. Schematic overview of a typical FCC unit is shown in Figure 1.
Ind. Eng. Chem. Res., Vol. 36, No. 2, 1997 391
Figure 2. MPC with a state estimation control loop.
Figure 1. Schematic overview of a FCC unit.
3. NLMPC with State Estimation In this section the structure of the MPC version developed by Ali and Zafiriou (1993) that utilizes directly the nonlinear model for output prediction is presented. Incorporating an optimization technique to estimate the modeling error on-line is also addressed. The objective of this paper is to implement a controller that can handle the fundamental control problem of FCC units rather than the testing of the performance of NLMPC through simulation. A usual MPC formation solves the following on-line optimization: P
min
|Γ(y(tk+i) - r(tk+i))|2 + ∑ ) i)1
∆u(tk),...,∆u(tk+M-1
M
|Λ∆u(tk+i-1)|2 ∑ i)1
(1)
filtering to perform on-line state estimation. The purpose is to reset observer states in the presence of modeling error or unmeasured disturbances. Detailed information on incorporating a Kalman filter with the NLMPC algorithm is given by Ali and Zafiriou (1993). However, the solution of the Lyapunov equation, which is required to compute the steady-state Kalman filter gain, requires the dynamic matrix of the linearized model to be stable (Lewis, 1986). Unfortunately, this is not the case in this specific example. In the present case, with the above characteristics, the state estimation technique is modified in a way that avoids the aforementioned problem. One way to deal with this situation is to introduce a completely arbitrary disturbance, say θ, in the right-hand side of eq 3 and then try to estimate it on-line. However, since θ is arbitrary, it is not easy to determine its end limits. Therefore, alternatively the unmeasured disturbances are estimated in a way which is analogous to the correction term computed by the Kalman filter (Ali and Zafiriou, 1993); i.e., the description of the unmeasured disturbance should be written as follows:
subject to
d ) θ[yˆ (tk) - g(x(tk))] AT∆U(tk) e b
(2)
For nonlinear MPC the predicted output y over the prediction horizon P is obtained by the numerical integration of:
x˘ ) f(x,u,t)
(5)
where the unknown parameter θ represents the steadystate Kalman gain for a continuous system and yˆ is the output measurement. θ, a diagonal matrix of size ny × ny, can then be determined on-line by solving the following objective function:
(3)
min |yˆ (tk) - y(tk)|2
(6)
0 e θi e 1, i ) 1, ..., ny
(7)
θ1,...,θny
y ) g(x)
(4) subject to
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 vector norm, and k denotes the current sampling point. Γ and Λ are diagonal weight matrices. r is the desired output trajectory. ∆U(tk) ) [∆u(tk) ... ∆u(tk+M-1)]T is a vector of M future changes of the manipulated variable vector u that are to be determined by the on-line optimization. A disturbance estimate should also be added to y in eq 1, or alternatively it can be absorbed in r. The latter is assumed in this work for the sake of simplicity. In the implementation in this present work, the disturbance is assumed constant over the prediction horizon and set equal to the difference between the plant and model outputs at present time k. The function of the “additive” constant disturbance in the model prediction is to introduce integral action and thus removes steadystate offset in the presence of model uncertainty or unmeasured disturbances. However, for open-loop unstable systems and/or slow dynamic processes, the model and the actual states may diverge away from each other and this thus degrades the regulatory performance of MPC. For this reason, the NLMPC formulation of Ali and Zafiriou (1993) also incorporates Kalman
where θi is the diagonal element of θ, ny is the number of measured outputs, and y is obtained by numerically integrating the state equations (3) and (4) augmented by the estimated disturbance. The augmented states are written as follows:
x˘ um ) fum(xms,xum,u)
(8)
x˘ ms ) fms(xms,xum,u) + d
(9)
x˘ ) [x˘ ms
x˘ um]T
y ) g(x)
(10) (11)
xms denotes the measured states and xum denotes the unmeasured states. Diagonal θ is used to reduce the number of parameters to be estimated. Although the trivial simplification is not optimal, it provides simple treatment against unmeasured disturbances and modeling errors. The mechanism of solving the MPC control problem is illustrated by Figure 2 and explained by the
392 Ind. Eng. Chem. Res., Vol. 36, No. 2, 1997
following algorithm: at any sampling time k, given the plant measurement yˆ (tk), the model output y(tk), and the past input and model states u(tk-1) and x(tk-1). Step 1: Obtain the optimum value of θ by solving the optimization problem (6) which involves the following steps: Step 1.1: For a given value of θ, which is determined by the optimization software, compute d from the available values of the measured and model outputs at k using eq 5. Step 1.2: Integrate eq 10 numerically from tk-1, using u(tk-1) and x(tk-1), up to tk to obtain a new value of x(tk). Step 1.3: Determine the model output y(tk) using eq 11. Step 1.4: If the objective function (6) satisfies the optimality tolerance set by the optimization software, go to step 2. Otherwise, go back to step 1.1. Step 2: Using the obtained optimum value of θ and the current value of d, solve the optimization problem (1) which involves the following steps: Step 2.1: For specific M future values of manipulated variables set by the optimization software, integrate the state equation (10) numerically over the prediction horizon P. Step 2.2: Obtain P future predicted outputs using eq 11. Correct the predictions by adding d to them. Step 2.3: If the objective function (1) does not satisfy the optimality tolerance set by the optimization software, then go to step 2.1. Otherwise, go to step 3. Step 3: Implement the obtained optimal manipulated variable for one sampling point. Set k ) k + 1 and go back to step 1. Correcting the predicted outputs by adding d in step 2.2 follows the internal model control approach to account for the unmeasured disturbances (Eaton and Rawlings, 1990). It should be noted that, in this present work, numerical solution of the differential equations is carried out using the DASSL software package (Petzold, 1983), and the optimization problems are solved using a DNCONF routine provided by the IMSL library. 4. Open-Loop Behavior The investigation of the uncontrolled (open-loop) FCC has shown that multiple steady states exist over a wide range of parameters. The bifurcation diagram of the open-loop case is shown in Figure 3. This open-loop analysis is based on the operating conditions given in Table 1, which are taken from an actual industrial unit as mentioned in section 1. There are three steady-state branches: a branch of low-temperature steady state with almost no reaction taking place, a high-temperature steady-state branch with low gasoline yield, and an intermediate unstable steady-state branch with a moderate temperature and high gasoline yield. The industrial operating condition corresponds to a steady state on the unstable intermediate branch of Figure 3 (point A on Figure 3 having yfa ) 0.872, x1 ) 1.5627, and x4 ) 0.380 88). It is clear that this operating point does not correspond to the maximum gasoline yield. (The maximum gasoline yield for this industrial unit corresponds to point B, having yfa ) 1.4145, x1 ) 1.230 37, and x4 ) 0.438 755.) The alteration of the operating conditions to operate at the point of maximum gasoline yield has been presented by Elnashaie et al. (1995). Operating at this point corresponds to a percentage increase in gasoline yield equal to about 15%. On the other hand, working around this point may not
Figure 3. Open-loop bifurcation diagram: (a) Dimensionless reactor dense phase temperature vs dimensionless air feed temperature to regenerator. (b) Dimensionless reactor dense phase gasoline concentration vs dimensionless air feed temperature to regenerator [solid line, stable steady states; dot lines, unstable DHP (degenerate Hopf point); A, industrial operating point; B, maximum gasoline concentration point].
be desirable because it is very close to the critical static turning point (which defines the boundary between the region of three steady states and five steady states) and also the degenerate Hopf bifurcation point (DHP). This DHP may give rise in its vicinity to homoclinical, torus, or period-doubling bifurcations which might lead to choas (Planeaux and Jensen, 1986). 5. Closed-Loop Behavior The steady-state analysis in the previous section indicates that optimal operation of the reactor occurs at the unstable steady-state point. This means that, in order to maximize the yield, the gasoline concentration should be maintained at that optimal operating point. Equivalently, since the gasoline concentration cannot be measured easily, maximum yield can also be achieved through maintaining the reactor temperature at its corresponding value at the desirable steady-state point. Thus, the control objective of this process is to stabilize and keep the reactor temperature at its desired set point, which corresponds to maximum gasoline yield. The simulation runs in the next section demonstrate that this goal can be attained by implementing NLMPC with state estimation during startup operation and when the process is under the influence of external disturbance. To show robustness of the NLMPC algorithm, the case when the model has parametric uncertainty has also been considered. The same control objective will also be examined by classical PI control for comparison purposes. A sampling time of 10 s will be used for all simulations in the following section. It should be noted that in all simulations both the plant and model are operated at the operating conditions given in Table 1, which are taken from an actual industrial unit as mentioned earlier. 5.1. Simulation Results. Perfect Model Case. In the following simulations a perfect model is assumed; i.e., there is no model-plant mismatch. Both classical
Ind. Eng. Chem. Res., Vol. 36, No. 2, 1997 393
Figure 5. Manipulated variable response for setpoint change using P control at various values of kc. For definitions of legends, see Figure 4.
Figure 4. Setpoint change response using P control at various values of kc. (a,b) Dashed line: kc ) 1.6. Dash-and-dotted line: kc ) 3. (c,d) Solid line: kc ) 5. Long dashed line: kc ) 10. Dashand-double-dotted line: kc ) 20. Table 2. Steady State state
initial value
new value
x1 x2 x3 x4 x5 x6
1.5 1.6 0.4 0.1 0.99 0.99
1.5627 2.008 0.38 0.38 0.99 0.99
PI and NLMPC schemes are examined, and their performance is compared. Servo Problem. The servo control objective is to start up the reactor dimensionless temperature (x1) from an initial value of 1.5 to a desirable value of 1.5627 using the air feed temperature (yfa) as the manipulated variable. The initial value of the manipulated variable is 0.87. Upper and lower bounds of magnitude of 0.57 and 1.87, respectively, were imposed on the manipulated variable. The corresponding values of the other states at the initial and final operating points are listed in Table 2. The performance of a proportional controller for this servo problem is shown in Figure 4 at various values for the controller gain (kc). The resulting behavior of the manipulated variable is shown in Figure 5 for the same values of kc. It is obvious from Figure 4 that at small values of kc (less than 5) the controller cannot keep the reactor temperature at its desired point. In fact, the temperature runs away to a higher operating point, diminishing the gasoline production. At larger values (5 and above) the controller becomes stable and maintains the temperature within the setpoint vicinity. However, temperature offset increases with the increase of the controller gain. Nevertheless, this offset has a marginal effect on the gasoline yield. The controller performance at kc ) 5, which corresponds to the best performance among those for the other values of kc, is replotted in Figure 6 for fare comparison with the NLMPC performance. It should be noted that a proportional controller solely is enough to handle the servo problem of a non-self-regulating process as it is known that the P controller is capable of removing offset for such processes. In fact, our investigation, which is not shown for brevity, revealed that adding integral action worsens the controller performance.
Figure 6. Setpoint change response. Solid line: using NLMPC with a perfect model, M ) 1, P ) 1, Λ ) [0], Γ ) [1]. Dashed line: using P control with kc ) 5.
The NLMPC performance for the same servo problem is illustrated by Figure 6. The values of the NLMPC parameters for this simulation are M ) 1, P ) 1, Λ ) [0], and Γ ) [1]. The figure also shows the time response of the regenerator temperature (x2), the gas oil concentration (x3), and the gasoline concentration (x4). The corresponding time response of the manipulated variable is shown in Figure 7. As indicated by Figure 6a, NLMPC managed successfully not only to bring the reactor temperature to the new desirable steady state but also to stabilize it around that point. Maintaining the reactor temperature around the desired value favored the gasoline production rate indicated by the increase of the gasoline concentration shown in Figure 6d. Specifically, x4 increased 280% over its initial value of 0.1 until it reached and settled around its desired steady-state value. Figure 7 shows the necessary variation in the manipulated variable required to start up the reactor temperature. According to the MPC tuning guidelines, large values of P and Λ are suggested for providing a stable feedback response. However, the small values used in the test for P and Λ were good enough, which are also favorable since, in that course, unnecessary increases in the computations of on-line optimization that NLMPC requires are avoided. Ex-
394 Ind. Eng. Chem. Res., Vol. 36, No. 2, 1997
Figure 7. Manipulated variable response for setpoint change. For definitions of legends, see Figure 6. Table 3. PI Controller Settings parameter
loop 1
loop 2
kc τI
400 000 -180 000
10 20
actly the same values of the NLMPC parameters are used in the next simulations and are found to give reasonable performance. As demonstrated by Figure 6, the temperature closedloop response obtained by the P controller has a higher overshoot and slower recovery to the setpoint than that obtained by NLMPC. However, the responses of the gasoline concentration in the two cases were somewhat comparable. The only drawback of the P control is that it requires careful detuning to obtain a reasonable performance, while a good performance was obtained by NLMPC even for the most aggressive case, i.e., zero weight on the manipulated variable and small P and M. Regulator Problem. In the following we demonstrate the effectiveness of NLMPC to remove the effect of unmeasured disturbance on the air feed temperature while the process is controlled and operating at the desired open-loop unstable steady state. Like in the servo case, the NLMPC performance will also be compared to that of a PI controller. The air feed temperature is used as the source of unmeasured disturbances just for demonstration purposes. Tow manipulated variables, namely, the air flow rate and the catalyst circulation rate, i.e., (Fa, Fc), will be used to regulate the process. The steady-state values for Fa and Fc are 10 670 and 88 605, respectively. Upper and lower bounds of (3000 and (20 000, as deviations from steady-state values, are imposed on Fa and Fc, respectively. The reason behind using two manipulated variables is to optimize their use since the controller may demand large changes of these variables due to the very small process gains relating the controlled output to these manipulated variables. First, we examine the effect of using conventional PI controls. Unlike the servo case, a proportional plus integral controller has to be used in order to eliminate steady-state error. A two-loop single PI control system configuration is tested. One loop pairs the controlled output with Fa and the other with Fc. Although this simultaneous pairing is not common in classical PI algorithms, it is done for comparison with NLMPC. The controller settings for each loop are obtained by the Ziegler and Nichols method (Ziegler and Nichols, 1942) and listed in Table 3. Figure 8 demonstrates the PI
Figure 8. Disturbance rejection response using PI control with kc1 ) 400 000, τI1 ) 10, kc2 ) -180 000, and τI2 ) 20.
Figure 9. Disturbance rejection response for the perfect model case using NLMPC without state estimation. M ) 1, P ) 1, Λ ) [0], Γ ) 1. Solid line: plant open-loop response. Dashed line: plant closed-loop response. Dash-and-double-dotted line: model closed-loop response.
closed-loop simulation for 20% step change in the air feed temperature. The figure shows the feedback response of the reactor temperature x1, the gasoline concentration x4, and the two manipulated variables Fa and Fc. It is shown that the Z-N settings give stable and satisfactory responses but at the expense of input saturation, which is depicted in Figure 8d, as expected from such a control algorithm. The simulation is replotted in Figure 11 for fare comparison with the NLMPC simulation. Open-loop and NLMPC closed-loop simulations for the same above load change in the nominal value of the air feed temperature are illustrated by Figure 9. The solid line corresponds to the open-loop case and the dashed line to the closed-loop case without using the proposed state estimation. The closed-loop response of the model without using the proposed state estimation is also shown in the same figure by the dash-and-double-dotted line. Figure 10 demonstrates the corresponding time response of the manipulated variables (Fa and Fc) for
Ind. Eng. Chem. Res., Vol. 36, No. 2, 1997 395
Figure 11. Disturbance rejection response for the perfect model case. Solid line: plant closed-loop response using NLMPC with state estimation. Dashed line: model closed-loop response using NLMPC with state estimation. Dash-and-double-dotted: plant closed-loop response using PI with the same controller parameters in Figure 8. Figure 10. Manipulated variable response for step change in disturbance load using NLMPC without state estimation for the perfect model case. M ) 1, P ) 1, Λ ) [0], Γ ) 1.
the closed-loop case. The NLMPC tuning parameters are the same as before. As expected, the open-loop simulation indicated, due to the internal instability of the process, substantial drift from the desired operating point. Specifically, a step change in the air feed temperature led the reactor temperature to move to another stable steady state associated with a dramatic drop in the gasoline concentration. Obviously, the controller helped in reducing the drift substantially, as illustrated by the dashed lines. However, the controller could not bring the process back to its desirable operating point due to the propagating mismatch between the plant and the model as shown by the difference between the dashed and the dash-and-double-dotted lines. It should be noted that tuning of the MPC parameters was tried to improve the closed-loop performance. However, this attempt was not successful because the poor MPC performance is essentially due to the propagating model-plant mismatch. For this reason it is necessary to reset the model states using the model proposed state estimation described in section 3. The improved performance of NLMPC by using the estimation technique is illustrated by Figures 11 and 12. The solid and the dashed lines of Figure 11 correspond to the NLMPC closed-loop response for the plant and the model, respectively, while the dash-and-double-dotted line corresponds to the PI closed-loop response using the Z-N settings. It is obvious from the figure that the proposed state estimation minimized substantially the mismatch especially for the controlled output (x1) in comparison to that of Figure 9. The reduction in the mismatch helped the NLMPC in perfectly rejecting the influence of the unmeasured disturbance and stabilizing the process. A 1.6% increase in the gasoline concentration over its desired steady-state value is associated with this simulation due to the 20% increase in the air feed temperature. Figure 12 demonstrates the necessary variation
Figure 12. Manipulated variable response for step change in disturbance load for the perfect model case. Solid line: using NLMPC with state estimation. Dash-and-double-dotted line: using PI control with the same controller parameters in Figure 8.
of the manipulated variables, which is determined by the controller, to reject the disturbance and to stabilize the system. Comparison of the two controllers through Figures 11 and 12 indicated their ability to suppress the disturbance influence. However, the improved PI performance is obtained via the Z-N tuning approach, which is somewhat computationally exhausting due to its trialand-error nature. On the other hand, a good NLMPC performance is obtained at minimum consumption of
396 Ind. Eng. Chem. Res., Vol. 36, No. 2, 1997
Figure 14. Manipulated variable response for setpoint change using NLMPC for the imperfect model case. M ) 1, P ) 1, Λ ) [0], Γ ) 1. Solid line: without state estimation Dashed line: with state estimation.
Figure 13. Setpoint change response using NLMPC for the imperfect model case. M ) 1, P ) 1, Λ ) [0], Γ ) 1. Solid line: without state estimation. Dashed line: with state estimation.
the manipulated variables, indicative of its optimal structure. Imperfect Model Case. In actual applications, the model is never perfect. Therefore, to show robustness of the NLMPC, we repeat all the above simulations when the model has uncertain parameters. Specifically, the dimensionless activation energies γ1, γ2, and γ3 of the model are assumed to be 20% larger than those of the plant. Due to this modeling error and given the fact that we are operating the process around an unstable operating point, poor performance of the NLMPC is expected. Hence, the proposed estimation technique must be used in all of the following simulations in order to achieve a good performance of the NLMPC. It should be noted that the same simulation and tuning parameters used before will also be used here. It should also be noted that detuning the MPC algorithm was also tried in the following simulations, but it was found useless for the same reason mentioned in the perfect model case. Furthermore, the traditional PI controls will not be tested in this case, as its performance will remain the same as before. This is because the PI controller is not a model-based algorithm; i.e., it is applied directly to the process. Thus, modeling error has no effect on its performance. Servo Problem. Feedback response of the process for the same previous servo problem is shown by Figures 13 and 14. Figure 13 illustrates the closed-loop response of the states, while Figure 14 illustrates the manipulated variable response. The solid line represents the performance of NLMPC without using the estimation technique. As illustrated by the solid lines, as the reactor temperature reaches the desired value which is unstable, the temperature surpasses it to a higher stable steady state, resulting in a sharp drop in the gasoline concentration. The growing model-plant mismatch, which is not shown in the figure for brevity, due to the parametric error, hindered NLMPC from successfully manipulating the air feed temperature in order to achieve the servo objective. On the other hand, when the proposed state estimation is used, the mismatch is reduced and consequently the NLMPC was able to force the reactor temperature to return to its setpoint as demonstrated by the dashed line. As in the
Figure 15. Disturbance rejection response using NLMPC for the imperfect model case. M ) 1, P ) 1, Λ ) [0], Γ ) 1. Solid line: open-loop response. Dashed line: closed-loop response without state estimation. Dash-and-double-dotted line: closed-loop response with state estimation.
perfect model case, a 280% increase in the gasoline concentration over its initial value is also observed in this simulation. This increase is due to the increase of the reactor temperature, which confirms the steadystate analysis discussed in section 4. Regulator Problem. Simulation of the process states and the manipulated variables for the regulator problem in the presence of parametric modeling error is illustrated by Figures 15 and 16, respectively. In reference to Figure 15, the solid line denotes the open-loop response, the dashed line denotes the closed-loop response without state estimation, and the dash-anddouble-dotted line denotes the closed-loop response with state estimation. As for Figure 16, the solid line denotes the closed-loop response using state estimation while the dashed line denotes the closed-loop response without state estimation. A poor regulatory performance of NLMPC, which is almost similar to that of the open loop, was observed, although additive disturbance estimates are added on the output prediction. As mentioned earlier, the latter could improve the NLMPC performance if the model-plant mismatch is constant over the prediction horizon (i.e., steplike). However, due to the internal instability and modeling error, the mismatch propagates with time, which is not shown in
Ind. Eng. Chem. Res., Vol. 36, No. 2, 1997 397
Moreover, NLMPC is shown to be more easily implemented than traditional PI, as the latter requires careful tuning. Nomenclature
Figure 16. Manipulated variable response for step change in disturbance load using NLMPC for the imperfect model case. M ) 1, P ) 1, Λ ) [0], Γ ) 1. Solid line: with state estimation. Dashed line: without state estimation.
the figures for brevity, detaining the NLMPC from achieving good performance. However, NLMPC can be improved if the proposed state estimation is used as shown by the dash-and-double-dotted curve of Figure 15. It is obvious that NLMPC rejected the effect of the disturbance; thus, the reactor temperature was brought back to and kept at its desirable setpoint. The latter resulted, as expected, in maintaining the gasoline concentration at maximum value. As in the perfect case, a 1.6% increase in the gasoline concentration over its desired steady-state value is also observed due to the 20% increase in yfa. 6. Conclusions A nonlinear model predictive control algorithm with modified state estimation and a conventional PI controller are tested and compared for the control of a type IV industrial fluid catalytic cracking (FCC) unit. The units are described by nonlinear models based on the two theories of fluidization for bubbling fluidized beds. The model has been previously verified against industrial units. The investigation showed that stabilization of the reactor temperature at the desired operating point can be achieved by manipulating the air feed temperature for the startup case and by the air feed flow rate and the catalyst circulation rate for the regulator case. Both cases demonstrated that a maximum gasoline yield can be obtained if the reactor temperature is kept within the desired value. It was also found that for imperfect model case good NLMPC regulator and servo performance can only be achieved if the model states are reset on-line. Comparison of PI with NLMPC revealed the existence of potential for an improved control performance through the use of NLMPC. This is due to the ability of NLMPC to utilize the manipulated variable in optimal fashion.
a1 ) normalized group, s-1, a1 ) FCCps/(AIRHRCpf) a2 ) normalized group, s-1, a2 ) FGCpl/(AIRHRCpf) a3 ) normalized group, s-1, a3 ) FCCps/(AIRHRCpa) A ) constant matrix for the linear constraints AIG, AIR ) area of bed within clouds for regenerator and reactor, respectively, m2 b ) lower and upper bounds of the linear constraints B, BD ) modified space velocity for reactor and regenerator, respectively, s-1 Ca1f, Ca2f ) concentration of gas oil and gasoline in input stream, respectively, kg/m3 CB ) coke burning rate, kg of coke burned/kg of solid s Cm ) total amount of coke deposits necessary for complete deactivation/total amount of catalyst Cpa, Cps ) specific heats of air and catalyst, respectively, kJ/(kg K) Cpl, Cpf ) specific heats of liquid and vaporized gas oil, respectively, kJ/(kg K) CRD ) dimensionless group connected with the coke balance in the reactor, s-1, CRD ) FCCm/(AIRHRCa1f) CGD ) dimensionless group connected with the coke balance in the regenerator, s-1, CGD ) FC/(AIGHG) d ) disturbance estimates EHR ) dimensionless dynamic parameters in heat balance equations for reactor EHG ) dimensionless dynamic parameters in heat balance equations for regenerator EMR ) dimensionless dynamic parameters in coke balance equations for reactor EMG ) dimensionless dynamic parameters in coke balance equations for regenerator Ei ) activation energy for reaction i, kJ/(kg mol) FA ) air flow rate, kg/s FC ) catalyst flow rate kg/s FG ) gas oil flow rate, kg/s HR, HG ) height of reactor and regenerator, respectively, m HLR, HLG ) heat loss in reactor and regenerator, respectively HRCB ) dimensionless rate of heat absorption due to endothermic cracking reaction in the reactor Hv ) normalized heat of vaporization rate for gas oil, s-1 M ) control horizon ny ) number of controlled outputs P ) output prediction horizon r ) setpoint RC ) normalized coke make rate, s-1 tk ) time at sampling instant k u ) manipulated variable vector x ) state vector xum ) unmeasured state vector xms ) measured state vector X1 ) dimensionless reactor dense phase temperature X2 ) dimensionless regenerator dense phase temperature X3 ) dimensionless gas oil concentration X4 ) dimensionless gasoline concentration X5 ) catalyst activity within reactor X6 ) catalyst activity within regenerator X3f ) dimensionless gas oil concentration in the feed X4f ) dimensionless gasoline concentration in the feed y ) model output yˆ ) process output Yfa ) dimensionless air feed temperature Yf ) dimensionless reactor feed temperature Yv ) dimensionless gas oil vaporization temperature
398 Ind. Eng. Chem. Res., Vol. 36, No. 2, 1997 Greek Letters R1, R2, R3 ) dimensionless preexponential factors for the cracking reaction, R1 ) K10Ca1f, R2 ) K20, R3 ) K30Ca1f Rs, Rd ) dynamic parameters for catalyst adsorption capacity for gas oil and gasoline, respectively βc ) exothermic factor for coke combustion reaction ) Cm(-∆H)/(TrefCpaFa) ∆u ) vector of manipulated variable change ∆U ) vector of M future manipulated variable changes ) bed voidage γ1, γ2, γ3 ) dimensionless activation energy for cracking reaction, γi ) Ei/RTref Γ ) diagonal weight matrix on the predicted error Λ ) diagonal weight matrix on the manipulated variable Fa, Fl, Ff ) density of air, liquid gas oil, and vaporized gas oil, respectively, kg/m3 θ ) diagonal matrix for the estimated disturbance
Literature Cited Ali, E.; Zafiriou, E. Optimization-Based tuning of Nonlinear Model Predictive Control With State Estimation. J. Process Control 1993, 3, 97-107. Arbel, A.; Huang, Z.; Rinard, I.; Shinner, R. Partial Control of FCC Units: Input Multiplicities and Control Structures. AIChE Annual Meeting, St. Louis, MO, 1993. Balchen, J. G.; Ljungquist, D.; Strand, S. State-Space Predictive Control. Chem. Eng. Sci. 1992, 47, 787-867. Cutler, C. R.; Gusciora, P.; McAmis, J.; Sorensen, R. Experiences Applying DMC on a Model IV FCC. AIChE Annual Meeting, Miami, FL, 1992. Eaton, J. W.; Rawlings, J. B. Feedback Control of Nonlinear Processes using On-line Optimization Techniques. Comput. Chem. Eng. 1990, 14, 469-479. Edwards, W. M.; Kim, H. N. Multiple Steady States in FCC Unit Operation. Chem. Eng. Sci. 1988, 37, 1611-1623. Elnashaie, S. S. E. H.; El-Hennawi, I. M. Multiplicity of the Steady States in Fluidized Bed Reactors, IV-Fluid Catalytic Cracking. Chem. Eng. Sci. 1979, 34, 1113-1121. Elnashaie, S. S. E. H.; Elshishini, S. S. Comparison between different Mathematical Models for the Simulation of Industrial Fluid Catalytic Cracking (FFC) Units. Math. Comput. Mode. 1993a, 18, 91-110. Elnashaie, S. S. E. H.; Elshishini, S. S. Digital Simulation of Industrial Fluid Catalytic Cracking Units-IV Dynamic Behavior. Chem. Eng. Sci. 1993b, 48, 567-583. Elnashaie, S. S. E. H.; Abasaeed, A. E.; Elshishini, S. S. Digital Simulation of Industrial Fluid Catalytic Cracking Units-V. Static and Dynamic Bifurcation. Chem. Eng. Sci. 1995, 50, 1635-1643. Elshishini, S. S.; Elnashaie, S. E. E. H. Digital Simulation of Industrial Fluid Catalytic Cracking Units-I, Bifurcation and its implications. Chem. Eng. Sci. 1990a, 45, 553-559. Elshishini, S. S.; Elnashaie, S. E. E. H. Digital Simulation of Industrial Fluid Catalytic Cracking Units-II, Effect of Charge Stock Composition on Bifurcation and Gasoline Yield. Chem. Eng. Sci. 1990b, 45, 2959-2964. Grosdidier, P.; Mason, A.; Aitolahi, A.; Heinonen, P.; Vanhamak, V. FCC Unit Reactor-Regenerator Control. Comput. Chem. Eng. 1993, 17, 165-179. Hidalgo, P. M.; Brosilow, C. B. Nonlinear Model Predictive Control of Styrene Polymerization at Unstable Operating Point. Comput. Chem. Eng. 1990, 14, 481-494. Huang, Z.; Spare, A. V.; Tsiligiannis, C. Dynamic modeling and control of fluidized catalytic cracking units. AIChE Annual Meeting, San Francisco, CA, 1989.
Iscol, L. The Dynamics and Stability of a Fluid Catalytic Cracker. Proceedings of the American Control Conference, Atlanta, GA, 1970; Paper 23-B. Jang, S.; Joseph, B.; Mukai, H. Control of Constrained Multivariable Nonlinear Process Using a Two-Phase Approach. Ind. Eng. Chem. Res. 1987, 26, 2106-2114. Kent, Z. Q.; Fisher, D. G. Model Predictive Control For Open-loop Unstable Process. Proceedings of the American Control Conference, San Francisco, CA, 1993; 791-795. Khandalekar, P. D.; Riggs, J. B. Nonlinear Process Model Based Control and Optimization of a Model IV FCC Unit. Comput. Chem. Eng. 1993, 19, 1153-1168. Kunii, D.; Levenspiel, O. Fluidization Engineering; John Wiley & Sons: New York, 1969. Kurihara, H. Optimal Control of Fluid Catalytic Cracking Processes. Ph.D. Dissertation, MIT, Cambridge, MA, 1967. Lee, W.; Kugelman, A. M. Number of Steady States Operating Points and Local Stability of Open-loop Fluid Catalytic Cracker. Ind. Eng. Chem. Process Des. Dev. 1973, 12, 197-204. Lee, J. H.; Morari, M.; Garcia, E. State Space Interpretation of Model Predictive Control. Automatica 1994, 30, 707-717. Lewis, F. L. Optimal Estimation with Introduction to Stochastic Control Theory; John Wiley & Sons: New York, 1986. Li, W. C.; Biegler, L. T. Newton-Type Controllers for Constrained Nonlinear Processes with Uncertainty. Ind. Eng. Chem. Res. 1990, 29, 1647-1657. Lokesh, K.; Georgakis, C. Effect of Process Nonlinearity on the Performance of Linear Model Predictive Controllers for the Environmentally Safe Operation of a Fluid Catalytic Cracking Unit. Ind. Eng. Chem. Res. 1994, 33, 3063-3069. Luyben, W. L.; Lamb, D. E. Feed-forward Control of a Fluidized Catalytic Reactor-Regenerator System. Chem. Eng. Prog. Symp. Ser. 1963, 59, 165-171. Petzold, L. R. A description of DASSL: A Differential-Algebraic System Solver. In Scientific Computing; Stepleman, R. S., Ed.; North-Holland: Amsterdam, The Netherlands, 1983. Planeaux, J. B.; Jensen, K. F. Bifurcation Phenomena in a CSTR Dynamics: a System with Extraneous Thermal Capacitance. Chem. Eng. Sci. 1986, 41, 1497-1523. Prett, D. M.; Gillette, R. D. Optimization and Constrained Multivariable Control of a Catalytic Cracking Unit. Proceedings of the American Control Conference, San Francisco, CA, 1980. Schmid, C.; Biegler, L. T. Application of Multistep Newton-Type Controllers To Fluid Catalytic Cracking. Proceeding of the American Control Conference, San Francisco, CA, 1990; pp 581-586. Weekman, V. W. A Model of Catalytic Cracking Conversion in Fixed, Moving and Fluid Bed Reactors. Ind. Eng. Chem. Process Des. Dev. 1968, 7, 90-97. Weekman, V. W.; Nace, D. M. Kinetics of Catalytic Cracking Selectivity in Fixed, Moving and Fluid Bed Reactors. AIChE 1970, 16, 397-404. Wright, G. T.; Breedjk, T.; Edgar, T. F. On-line Parameter Estimation and Adaptation in Nonlinear Model-Based Control. Proceedings of the American Control Conference, Boston, 1991. Ziegler, J. G.; Nichols, N. B. Optimum Settings of Automatic Controllers. Trans. ASME 1942, 64, 759.
Received for review June 21, 1996 Revised manuscript received October 18, 1996 Accepted October 21, 1996X IE9603575
X Abstract published in Advance ACS Abstracts, December 1, 1996.