- Email: [email protected]
Nonlinear model predictive control of a reactive distillation column
ARTICLE IN PRESS
Control Engineering Practice 15 (2007) 231–239 www.elsevier.com/locate/conengprac
Nonlinear model predictive control of a reactive distillation column Rohit Kawathekar, James B. Riggs Chemical Engineering, Texas Tech University, Lubbock, TX 79409-3121, USA Received 5 July 2005; accepted 10 July 2006 Available online 7 September 2006
Abstract This paper considers the application of nonlinear model predictive control (NLMPC) to a highly nonlinear reactive distillation column. NLMPC was applied as a nonlinear programming problem using orthogonal collocation on finite elements to approximate the ODEs that constitute the model equations for the reactive distillation column. Diagonal PI controls were used to identify that the [L/D,V] and the [L/D,V/B] configurations performed best. NLMPC was applied using the [L/D,V] configuration and found to provide a factor of 2–3 better performance than the corresponding PI controller. The effect of process/model mismatch on the performance of the NLMPC controller was also evaluated. r 2006 Elsevier Ltd. All rights reserved. Keywords: Reactive distillation; Nonlinear model predictive control; Configuration selection; Process/model mismatch
1. Introduction Reactive distillation combines both separation and reaction in one unit and has been applied industrially for a number of years. Reactive distillation can offer significant economic advantages for certain cases, particularly for systems that involve reversible reactions. Examples of commercial successes of reactive distillation include nylon 6,6 processes, methyl acetate processes, ethyl acetate processes and methyl tert-butyl ether processes (Doherty & Buzad, 1992). Most of the literature available on reactive distillation is based on steady-state conditions including process design (e.g., Chang & Seader, 1988; Krishnamurthy & Taylor, 1985) and the analysis of multiple steady states (e.g., Abufares & Douglas, 1995; Jacobs & Krishna, 1993; Sneesby, Tade, & Smith, 1998). Dynamic modeling and simulation (e.g., Grosser, Doherty, & Malone, 1987; Ruiz, Basualdo, & Scenna, 1995) have also been studied, but a relatively small amount of work has been reported on the control of reactive distillation columns. Chen, Hontoir, Huang, Zhang, and Julian Morris (2004) performed a nonlinear dynamic modeling case study of a reactive Corresponding author. Tel.: +1 806 742 1765; fax: +1 806 742 3552.
E-mail address: [email protected] (J.B. Riggs). 0967-0661/$ - see front matter r 2006 Elsevier Ltd. All rights reserved. doi:10.1016/j.conengprac.2006.07.004
distillation column using a nonlinear first-principles model combined with an artificial neural network model. Reactive distillation is a challenge for control due to process nonlinearity and complex interactions between vapor–liquid equilibrium and chemical reactions. Al-Arfaj and Luyben (2000), Sneesby, Tade, Datta, and Smith (1997) and Kumar and Daoutidis (1999) have considered decentralized PI control structures for reactive distillation columns. Khaledi and Young (2005) applied PI controls and linear model predictive control to a dynamic simulation of an ETBE reactive distillation column using the (L,V) configuration. The linear model predictive controller used first-order plus deadtime models for each of its input/output models. Kumar and Daoutidis (1999) have discussed the superior performance of nonlinear controllers compared to linear controllers for reactive distillation systems. Vo¨lker, Sonntag, and Engell (2006) applied a linear reduced-order MPC controller to a semibatch reactive distillation column. Gru¨ner et al. (2003) developed a general controller for reactive distillation systems based on asymptotically exact input/output linearization and applied their controller to a simulation of an industrial reactive distillation column with improved performance compared to a well-tuned linear controller. This paper considers the application of nonlinear model predictive control (NLMPC) to an ethyl acetate reactive
ARTICLE IN PRESS R. Kawathekar, J.B. Riggs / Control Engineering Practice 15 (2007) 231–239
232
Nomenclature
xi,j
B D Fi Hi hi ki
yi,j
kfi L L/D Mi l nc Pi Psi,j QC QR QMS ri u V V/B Voli
bottom product flow rate distillate product flow rate the external feed rate to the ith stage. the enthalpy of the vapor leaving the ith stage the enthalpy of the liquid leaving the ith stage the reaction equilibrium rate constant for the ith stage the reaction rate constant for the ith stage the reflux rate the reflux ratio the molar holdup of liquid on the ith stage the number of components the pressure on the ith stage the vapor pressure of component j on the ith stage the condenser duty the reboiler duty the mover suppression factor for the NLMPC controller the reaction rate on the ith stage the manipulate variable the boilup rate the boilup ratio the volume of the liquid holdup on the ith stage
distillation column. In Section 2, the case study is described in detail. Application of diagonal PI controls is presented in Section 3. Development and implementation of NLMPC algorithm is discussed in Section 4. In Section 5, the comparison between PI and NLMPC control results is presented. The effect of process-model mismatch on the closed-loop performance of NLMPC controller is analyzed in Section 6. 2. Model development This study is based on the production of ethyl acetate by a reactive distillation using the esterification reaction between acetic acid and ethanol. The achievable conversion in this reversible reaction is limited by the equilibrium conversion. The reaction is slightly endothermic and takes place in the liquid phase. Though the esterification reaction is self-catalyzed, sulfuric acid can act as external catalyst to enhance the reaction rate: CH3 COOH þ C2 H5 OH Sulfuric acid
!
CH3 COOC2 H5 þ H2 O:
ð1Þ
The reaction kinetics for this reaction were proposed by Aljeski and Duprat (1996) and include reaction equilibrium limitations. Bock and Wozny (1997) reported a detailed analysis of ethyl acetate reactive distillation column. The slow reaction rate for this system leads to lower conversion
ysp
the mole fraction of the jth component in the liquid on the ith stage the mole fraction of the jth component in the vapor on the ith stage the setpoint for the controlled variable
Greek Symbols gi,j mj ri
the liquid phase activity coefficient for the jth component on the ith stage the stochiometric coefficient for the jth component in the primary reaction the liquid density for the ith stage
Subscripts A B C c D i f j n
acetic acid ethyl alcohol ethyl acetate condenser water stage number feed component number the top tray number
for the reactive distillation column. A reactive column with countercurrent flow of reactants was used. Pure ethanol and acetic acid were fed separately into a column that operates at atmospheric pressure. Under these conditions, acetic acid is the heaviest of the components and moves toward the bottom of the column. Ethyl acetate is the lightest and moves toward the top of the column. It is expected that the middle portion of the distillation column is the primary reaction zone. The rectifying section fractionates the ethyl acetate out of acetic acid, and the stripping section removes alcohol from water. This quaternary system consisting of ethanol, acetic acid, water and ethyl acetate is highly nonideal. It can form four binary azeotrope mixtures and one ternary azeotrope. Over a wide range of composition, ethanol and water do not differ greatly in volatility, making it difficult to produce only water as the bottom product. The reactant ethanol has a relatively high volatility in the reaction zone. This leads to a low composition of ethanol in the liquid phase, reducing the production rate of ethyl acetate. Under these conditions, the rate of the esterification reaction between acetic acid and ethanol is generally low, which implies that it is favored by long residence times in each stage. It is evident from all the previous studies on ethyl acetate reactive distillation columns that an unfavorable physical equilibrium makes the production of highpurity ethyl acetate impossible from a single distillation column.
ARTICLE IN PRESS R. Kawathekar, J.B. Riggs / Control Engineering Practice 15 (2007) 231–239 Reactive Column
AA feed
Recovery Column
101.3 kPa
350.3 kPa Purge
EtOH feed
RX
RC
Ethyl acetate product Recycle
Fig. 1. Two-column sequence for the production of high-purity ethyl acetate.
However, Seferlis and Grievink (2001) proposed adding a recovery column (Fig. 1), which operates at a higher pressure to produce a high-purity ethyl acetate product. The recovery column was combined with the reactive distillation column and designed to operate at a higher pressure (350 kPa) to break the azeotropes, which increases the overall conversion and produces a high-purity ethyl acetate product. The distillate stream from the reactive column is fed to the recovery column. At a relatively high pressure, ethyl acetate becomes heavier than ethanol and water so that it appears as the bottom product. The target purity level of the ethyl acetate was set at 99.5%. Reaction in the recovery column is negligible because the column operates without any sulfuric acid catalyst and the liquid holdup on the trays is small. Virtually all of the acetic acid in the overhead of the reactive column ends up in the bottom product of the recovery column, and therefore, directly affects the purity of final ethyl acetate product. Hence, a specification is imposed on the maximum allowable concentration of acetic acid in the distillate of the reactive column. To assist in the reduction of the acetic acid in the distillate of the reactive column, the reactive column was designed with a reduced liquid phase holdup in the upper section of trays of the reactive column to suppress the reverse reaction of ethyl acetate. Seferlis and Grievink (2001) developed an economically optimized steady-state design for a reactive distillation process based on the two-column configuration (Fig. 1). The exact design and steady-state operating conditions of the reactive distillation column from Seferlis and Grievink (2001) were used as the basis of this study while the recovery column was not directly modeled. An equilibrium stage dynamic model of the reactive distillation column was developed and benchmarked against the results
233
reported by Seferlis and Grievink (2001). This dynamic model served as the process that was controlled by the PI and NLMPC controllers in this study. An equilibrium stage dynamic model of the reactive distillation column was implemented with the following assumptions: (1) each tray was assumed to be an ideal stage, (2) nonequal molar overflow was considered, (3) the liquid activity coefficient was modeled using the approach of Suzuki, Yagi, Komatsu, and Hirata (1971) [15], (4) the SRK method (Soave, 1992) was used to model the enthalpy departure functions, (5) a hydraulic time constant was used to model the liquid flow from each tray, (6) tray temperatures in the rectifying and stripping sections were used to infer the impurity composition of the overhead and bottom products, respectively, (7) the heat of dilution of acetic acid in water was neglected and (8) equilibrium limited reaction kinetics (Aljeski & Duprat, 1996) are assumed. For simulation purposes, the process is modeled using a detailed tray-to-tray model. For a standard reactive distillation tray where vapor holdup is considered negligible and reaction takes place in liquid phase, the modeling equations can be written as follows: nc X dðM li Þ ¼ F i þ Liþ1 þ V i1 Li V i þ Vol i m j ri , dt j¼1
dðM li xi;j Þ ¼ F i xi;f þ Liþ1 xiþ1;j þ V i1 yi1;j dt Li xi;j V i yi;j þ Vol i mj ri ,
(2)
ð3Þ
dðM li hi Þ ¼ F i hf i þ Liþ1 hiþ1 V i1 H i1 Li hi V i H i , dt
(4)
dhi 0, dt
(5)
1 ri ¼ r2i kfi xA;i xB;i xC;i xD;i , ki
(6)
yi;j Pi ¼ gi;j xi;j Psi;j .
(7)
For the partial condenser, the material balance and energy balance equations are as follows: nc X dðM lc Þ ¼ V n Lc V D þ Vol c m j rc , dt j¼1
(8)
dðM lc xc;j Þ ¼ V n yn;j Lc xc;j V D yD;j þ Vol c mj rc , dt
(9)
dðM lc hc Þ ¼ V n H n Lc hc V D H D Qc , dt
(10)
dhlc 0, dt
(11)
ARTICLE IN PRESS R. Kawathekar, J.B. Riggs / Control Engineering Practice 15 (2007) 231–239
234
ri ¼
r2i kfi
1 xA;i xB;i xC;i xD;i , ki
(12)
yi;j Pi ¼ gi;j xi;j Psi;j ,
(13)
where subscript c represents condenser stage, subscript n represents the top tray of the column and subscript D represents the distillate. VD represents the vapor distillate from the partial condenser stage. The reflux accumulator stage is modeled as a liquid holdup tank. The terms relating to the vapor phase, i.e., vapor flow, vapor composition, vapor enthalpy corresponding to reflux accumulator are absent in the material and energy balances. The material and energy balance equations for the reboiler are written as follows: nc X dðM lB Þ ¼ Li B V B þ Vol B mj rB , dt j¼1
(14)
dðM lB xB;j Þ ¼ Li xi;j BxB;j V B yB;j þ Vol B mj rB , dt
(15)
dðM lB hB Þ ¼ Li hi BhB V B H B þ QR , dt
(16)
dhlB 0, dt
(17)
1 ri ¼ r2i kfi xA;i xB;i xC;i xD;i , ki
(18)
yi;j Pi ¼ gi;j xi;j Psi;j .
(19)
D
PC
LT
RSP LC
CW LC
LT PT
L
Fig. 2. Schematic of the column pressure controls and the accumulator level controls.
The overhead level and pressure control structures are shown in Fig. 2 based on producing a vapor distillate product. The overhead pressure control loop sets the distillate vapor flow rate. The overhead pressure was modeled considering the total number of moles of vapor in the overhead and the vapor volume in the overhead line and the overhead condenser. A flooded condenser is used and the level of liquid in the condenser is used to vary the surface area in the condenser available for condensing the overhead vapor. An accumulator level controller determines the setpoint for the condenser level controller. Because all the acetic acid that leaves in the overhead of the reactive distillation column ends up in the ethyl acetate product produced by the recovery column, the primary composition control point is the acetic acid concentration in the overhead product from the reactive distillation column (i.e., the design setpoint is 0.04% acetic acid in the overhead). The composition setpoint for the bottoms is 0.08% ethyl acetate, which was determined by Seferlis and Grievink (2001) using a steady-state economic analysis of the two-column system. The temperature of the second tray from the top was used to infer the impurity level in the overhead product and the fifth tray from the bottom was used to infer the impurity level in the bottom product. These tray temperatures were identified using a steady-state model of the reactive distillation column to determine the tray temperatures that exhibited the strongest correlation between changes in the product purities and changes in the tray temperatures. The correlations between tray temperature and product composition were updated using new composition analyzer information when it became available. It was determined that this inferential tray temperature scheme provided excellent estimates of the product composition by comparing the inferential estimate of the product impurity levels with the corresponding values generated by the process simulator. Therefore, an observer was not considered necessary in this case. Table 1 lists the steady-state design conditions for the reactive distillation column studied here. The model equations for the reactive distillation column (i.e., the liquid dynamics equation, the component material balances and the energy balance for each tray, the component material balances and the overall material balances for the accumulator, reboiler and condenser, and the PI control equations for the pressure controller and the level controllers) were integrated numerically using the Euler method with a step size of 0.5 s. To improve the computational efficiency of the dynamic reactive distillation model, the inside-out algorithm (Boston & Sullivan, 1974) was used to calculate the component K-values. If a tray temperature changed by more than 1 1C, the parameters of the K-value correlation were updated. In addition, all K-value correlations were updated each 30 s of simulated time. The accumulator and reboiler level controllers were tuned for critically damped behavior.
ARTICLE IN PRESS R. Kawathekar, J.B. Riggs / Control Engineering Practice 15 (2007) 231–239
4.6E-04
31 25–3a 7–24 1–6
Stage holdups Rectifying section (m3) Reactive section (m3) Stripping section (m3)
0.151 1.75 1.65
Distillate flow (kmol/s) Bottoms flow (kmol/s) Reflux ratio
4.88 103 9.09 104 0.8211
Feed flow rate Acetic acid feed (kmol/s) Ethanol feed (kmol/s) Recycle feed (kmol/s)
8.92 104 8.74 104 4.03 103
Distillate composition (mol frac) Ethanol Acetic acid Water Ethyl acetate
0.2383 0.0004 0.2093 0.5520 0.0249 0.0876 0.8867 0.0008 171.0
Condenser duty (kW)
137.7
Low holdup was used to reduce reaction in the rectifying section.
Table 2 Steady-state changes for a 10% change in a manipulated variable Configuration
[L/D,V] [L/D,V/B] [L,V] [L,V/B]
L/D,V L/D,V/B L/D,B
4.4E-04 4.2E-04 4.0E-04 3.8E-04 0
500
1000
1500
2000
Time (min) Fig. 3. Overhead composition control using PI controls for an unmeasured feed composition upset.
3. Diagonal PI control results
Bottoms composition (mol frac) Ethanol Acetic acid Water Ethyl acetate Reboiler duty (kW)
a
Overhead impurity (mole fraction)
Table 1 Steady-state design conditions for the reactive distillation column Number of stages Rectifying section Reactive section Stripping section
235
% change in the process gain Overhead
Bottom
200 150 150 33
300 300 600 300
The steady-state gains for the ethyl acetate reactive distillation column are listed in Table 2. Note that for approximately a 10% change in each manipulated variable of each configuration, the process gain mostly changed by 150–600%, indicating severe steady-state process nonlinearity. In addition, several of the off-diagonal gains demonstrated bi-directionality. As is usually the case, process integration (i.e., adding the recovery column with recycle in this case) offers significant steady-state economic benefits but at the expense of a much more challenging control problem. The dynamic model of the reactive distillation column described in this section was used as the process, which was controlled by PI and nonlinear MPC controllers.
Diagonal PI composition controllers were used to compare the control performance of various control configurations and as a benchmark for the NLMPC control results. Because the primary control objective for the reactive distillation column is control of the impurity in the overhead product, the bottom composition control loop was tuned first for a closed-loop overdamped response. Then the overhead PI composition controller was tuned using setpoint changes. The overhead tuning procedure was based on using an ATV test using a single tuning factor (Riggs, 2001) combined with Tyreus– Luyben settings (Tyreus & Luyben, 1992) and finally adjusted (Riggs, 2001) to minimize the IAE from setpoint for a series of setpoint changes. The control performance of the overhead and bottoms composition controllers was evaluated for several unmeasured disturbances in feed rate and feed composition. Because fresh acetic acid and ethyl alcohol are each fed to the reactive column as well as the recycle stream from the recovery column, feed composition upsets were implemented by changing the ratio of fresh acetic acid and ethyl alcohol and feed rate changes were implemented by changing the recycle feed rate. The [L,B], [L,V], [L,V/B], [L/D,B], [L/D,V] and [L/D,V/ B] configurations were compared for control performance based on dual composition PI control (note that for a [A,B] configuration A is manipulated to control the overhead composition and B to control the bottom composition). Configurations using the distillate rate, D, as a manipulated variable were not considered because D is used to control the column pressure (Fig. 2). Figs. 3 and 4 show the overhead and bottoms composition control performance for the [L/D,B], [L/D,V] and [L/ D,V/B] configurations for an unmeasured feed composition upset. Table 3 lists the IAEs for each configuration for a feed rate upset. Based on these results, the [L/D,V] and [L/ D,V/B] configurations provided the best overall control performance.
ARTICLE IN PRESS R. Kawathekar, J.B. Riggs / Control Engineering Practice 15 (2007) 231–239
236
differences observed between the process output and the model prediction are due to additive step disturbances in the output. These disturbance terms are also assumed to remain constant over the prediction horizon.
Bottom ethyl acetate impurity (mole fraction)
8.3E-04 8.2E-04 8.1E-04 8.0E-04 7.9E-04
4.1. Solution algorithm
7.8E-04 7.7E-04
L/D,V/B
7.6E-04
L/D,V
7.5E-04
L/D,B
7.4E-04 0
2000
4000
6000
Time (min)
Fig. 4. Bottom composition control performance using PI control for an unmeasured feed composition upset.
Table 3 Dual PI composition control performance indices for unmeasured feed rate disturbance Configuration
IAE
[L,V] [L,B] [L,V/B] [L/D,V] [L/D,B] [L/D,V/B]
Overhead
Bottoms
0.143 0.177 0.142 0.115 0.393 0.114
1.278 0.830 1.239 1.270 1.215 0.988
4. Implementation of NLMPC NLMPC uses nonlinear models in a model predictive control framework to choose control action. The NLMPC controller chooses the future control moves that minimize the following objective function: min
N 1 X j¼0
½ysp yj 2 þ
N 1 X
QMS ½Duj 2 ,
(20)
j¼0
subject to bounds on manipulated, output and state variables. The first term represents the error from setpoint and the second term represents the weighted variation in the manipulated variables. Once the optimum set of future moves is determined, the initial control move is implemented. In the next control cycle, the next set of optimum control moves is determined. Even though the full set of optimum control moves into the future are calculated at each control interval, only the first move is actually implemented at each control interval. The formation of the NLMPC objective function will involve a disturbance term, which is added to the output prediction over the entire prediction horizon to match the model prediction to the current process measurement. The NLMPC procedure in this work assumes that the
NLMPC was implemented by setting up the control problem as a nonlinear programming (NLP) problem and solving it over the prediction horizon. It is necessary to simultaneously solve an optimization problem (the tradeoff between control to setpoint and changes in the manipulated variables) and the system model equations. These two procedures may be implemented either sequentially or simultaneously. 4.1.1. Sequential solution and optimization algorithm The sequential algorithm employs separate algorithms to solve the differential equations and to carry out the optimization. First, a manipulated variable profile is selected and the differential equations are solved numerically to obtain the controlled variable profile. The objective function is then determined. The gradient of the objective function with respect to the manipulated variables is determined either by numerical perturbation or by using analytical derivatives. The control profile is then updated using an optimization algorithm. The process is repeated until the optimal profile is obtained. This is referred as a sequential solution and optimization algorithm. The implementation of the sequential solution and optimization algorithm is relatively simple to apply. However, there are some drawbacks associated with this approach. The sequential solution and optimization requires the solution of differential equations at each iteration of the optimization. Jones and Finch (1984) found that such methods spend about 85% of the time integrating the model equations in order to obtain gradient information. This can make the implementation of this algorithm computationally expensive for cases involving a large number of model equations. The gradient information required for the optimization procedure is often obtained through numerical differentiation, as the analytical derivatives are not always available for highly nonlinear model equations involving complicated thermodynamic relations. To obtain the gradients using finite difference typically involves differencing the output of an integration routine with adaptive step sizes. Gill, Murray, and Saunders (1988) points out that the integration error is unpredictable and hence differencing output of an integration routine greatly degrades the quality of the finite difference derivatives. It is also difficult to incorporate the constraints on state variables with the use of the sequential solution and optimization approach (Meadows & Rawlings, 1997). 4.1.2. Simultaneous solution and optimization algorithm The simultaneous solution and optimization algorithm involves the model equations appended to the optimization
ARTICLE IN PRESS R. Kawathekar, J.B. Riggs / Control Engineering Practice 15 (2007) 231–239
problem as equality constraints. Then the NLP problem is posed to optimize the objective function such that the algebraic form of the model differential equations are satisfied, and the upper and lower constraints on states, controlled variables and manipulated variables are met. This can greatly increase the size of optimization problem, leading to a tradeoff between the two approaches. Meadows and Rawlings (1997) reported that for small problems with few states and a short prediction horizon, the sequential solution and optimization algorithm is more computationally efficient. For larger problems, the simultaneous solution and optimization approach is more computationally efficient and more reliable. In the simultaneous solution and optimization approach, the model equations and the process constraints (e.g., a manipulated variable cannot be less than zero) become constraints for the optimization problem to determine the future control moves. Since the model equations largely appear as ODEs, orthogonal collocation on finite elements (Finlayson, 1980) was used to convert the ODEs in time into a set of algebraic equality constraints for the optimization problem. Due to the size of the resulting NLMPC problem considered here, the simultaneous solution approach was chosen for this work. The process model and the NLMPC controller were implemented in Fortran code. The NLP problem resulting from the implementation of NLMPC was solved using the SNOPT software (Gill et al., 1988). Typical computational times for 30 h of simulation of the NLMPC controller using an 1100 MHz in an AMD PC with Windows XP required 10 h CPU. 4.2. Selection of tuning parameters The tuning parameters that have a significant effect on NLMPC performance are the prediction horizon, control horizon, sampling interval, penalty weight matrices and move suppression factors. A set of heuristics based on the linear systems and numerical simulations were used to select the final tuning parameters Following is a discussion of the tuning parameters used for the tuning of each NLMPC controller considered: (1) The control interval was set at 20 min, which was shown to be equivalent to continuous control using a correlation by Marlin (1995). (2) The prediction horizon was set at 80 by using different levels and comparing control performance. (3) A control horizon of 15 was found to provide good control performance. (4) The equal concern errors were selected to provide a 10 times greater weighting of the overhead errors from setpoint than the bottoms because impurities in the overhead directly affect the purity of the ethyl acetate final product. (5) The move suppression factors were adjusted to minimize the IAE from setpoint for the overhead product for a series of overhead impurity setpoint changes. 5. Comparison between PI and NLMPC The NLMPC control algorithm described in the previous sections was applied for dual-ended composition control of
237
Table 4 Ethyl acetate NLMPC control performance indices for overhead impurity setpoint tracking Configuration
Overhead loop IAE
Bottoms loop IAE
[L/D,V] NLMPC [L/D,V] PI
0.91 2.41
2.14 6.84
Table 5 NLMPC control performance indices for unmeasured feed composition disturbance rejection Configuration
[L/D,V] NLMPC [L/D,V] PI
IAE Overhead
Bottoms
0.0609 0.115
0.483 1.270
the simulation of the ethyl acetate reactive distillation column. The previous PI dual-ended composition control results indicated that [L/D, V] configuration provided good control performance for setpoint tracking as well as the rejection of unmeasured disturbances. NLMPC was applied using the [L/D, V] configuration to determine the benefits compared to conventional PI controls. Table 4 summarizes the comparisons between diagonal PI control and NLMPC using the [L/D,V] configuration for a series of setpoint changes in the overhead impurity. Table 5 summarizes the comparisons between diagonal PI control and NLMPC for an unmeasured feed composition upset. These results are based on an NLMPC controller with an almost perfect process model. Based on the IAEs of these tests, the NLMPC controller outperformed the PI controller by a factor of 2–3. Due to the severe nonlinearity of this process, one would expect a larger disparity between PI and NLMPC. On the other hand, this column was operated largely as a single-ended column control problem because the bottom composition control was detuned compared to the more important overhead product, which significantly simplified the control problem. 6. Process/model mismatch Process/model mismatch for the NLMPC controller was represented by introducing a 5% and a 25% error between the reaction equilibrium constant used in the dynamic process simulator and the model equations used by the NLMPC controller to calculate control action. Figs. 5–9 show the effect of process model/mismatch for NLMPC for an unmeasured feed composition upset. The PI control results are also shown for this case. To more completely analyze the sensitivity of PI and NLMPC control performance to unmeasured disturbances, a sinusoidal disturbance in the recycle feed rate with constant amplitude and varying frequency was applied to the process with a PI controller and an NLMPC controller
ARTICLE IN PRESS R. Kawathekar, J.B. Riggs / Control Engineering Practice 15 (2007) 231–239
238
1.68E+02 Reboiler duty (btu/sec)
Overhead impurity (mole fraction)
4.2E-04 4.0E-04 3.8E-04 0% mismatch 5% mismatch 25% mismatch PI
3.6E-04 3.4E-04 0
200
400
600
Time (min)
25% mismatch PI
1.62E+02 1.60E+02 1.58E+02 1.56E+02 1.54E+02 500
1500
0.07 PI
Overhead impurity amplitude ratio
0.06
8.1E-04
7.9E-04 0% mismatch 5% mismatch 25% mismatch PI
7.7E-04
NLMPC
0.05
NLMPC
0.04
5% mismatch
0.03 0.02 0.01 0 0
7.5E-04 0
500
1000
Time (min)
Fig. 6. Control results for the bottom product for an unmeasured feed composition upset showing the effect of process/model mismatch for NLMPC and compared with PI control for the [L/D,V] configuration.
1.3 1.2 1.1 1.0 0% mismatch 0.9
5% mismatch 25% mismatch
0.8
PI 0.7 0
1000 Time (min)
Fig. 8. The manipulated variable profile for bottoms for an unmeasured feed composition upset for the [L/D,V] configuration.
8.3E-04 Bottoms ethyl acetate impurity (mole fraction)
5% mismatch
1.64E+02
0
Fig. 5. Control results for the overhead product for an unmeasured feed composition upset showing the effect of process/model mismatch for NLMPC and compared with PI control for the [L/D,V] configuration.
Reflux ratio
0% mismatch
1.66E+02
200
400
600
Time (min)
Fig. 7. The manipulated variable profile for overhead for an unmeasured feed composition upset for the [L/D,V] configuration.
in use. A sinusoidal disturbance to the recycle feed rate was applied with an amplitude of 10% of the base case recycle feed rate. The closed-loop amplitude ratio for the overhead impurity was calculated by taking the ratio of the amplitude of overhead impurity resulting from the
500
1000
1500
Time period for sinusoidal feed disturbance (min)
Fig. 9. A comparison of closed-loop amplitude ratios for PI and NLMPC for the [L/D,V] configuration.
sinusoidal disturbance to the amplitude of the recycle feed disturbance. The closed-loop performances of the PI control structure, NLMPC controller as well as NLMPC with 5% process/model mismatch were compared using the unmeasured sinusoidal recycle feed rate disturbance. A plot of the closed-loop amplitude ratio versus the period of the sinusoidal feed disturbance (Fig. 9) shows that for the PI controller as well as the NLMPC controller exhibit sensitivity to the feed disturbance in a specific range of frequencies (120 min period to a 720 min period). The closed-loop amplitude ratio for the diagonal PI control structure is shown to be 2–3 times larger than that for NLMPC over the sensitive range. At low periods (high frequencies), the inertia of the process damps out the effect of the sinusoidal disturbance. At large periods (low frequencies), the feedback controllers have ample time to absorb the disturbances by feedback. Fig. 9 also indicates that the 5% process/model mismatch does not significantly affect the control performance for the NLMPC controller. These results are consistent with the previous results for the comparisons between PI and NLMPC control performance for the setpoint tracking as well as unmeasured disturbance rejection, but more quantitatively define the
ARTICLE IN PRESS R. Kawathekar, J.B. Riggs / Control Engineering Practice 15 (2007) 231–239
control performance differences between PI and NLMPC for this case. 7. Conclusions As expected, NLMPC was found to provide significantly better control performance than PI controls, but this work quantifies the advantage as a factor of 2–3 reduction in variability. The advantage of the NLMPC controller comes from faster closed-loop dynamic performance compared to the PI controller resulting from using a nonlinear dynamic model of the reactive distillation column considered here. In addition, NLMPC was found not to be particularly sensitive to process/model mismatch. Even though the reactive distillation column considered in this study was shown to be extremely nonlinear, the PI control was able to control the process reasonably well. Acknowledgements This work was supported by the member companies of the Texas Tech Process Control and Optimization Consortium. The authors would also like to acknowledge assistance from Charles R. Cutler, Cutler-Technologies. References Abufares, A. A., & Douglas, P. L. (1995). Mathematical modeling and simulation of an MTBE catalytic distillation process using SPEEDUP and aspen plus. Chemical Engineering Research and Design, Transactions of the Institute of Chemical Engineers Part A, 73, 3–12. Al-Arfaj, M. A., & Luyben, W. L. (2000). Comparative control of ideal and methyl acetate reactive distillation. Chemical Engineering Science, 57, 5039–5050. Aljeski, K., & Duprat, F. (1996). Dynamic simulation of the multicomponent reactive distillation. Chemical Engineering Science, 51, 4237–4252. Bock, H., & Wozny, G. (1997). Analysis of distillation and reaction rate in reactive distillation. Distillation and absorption ‘97. Institute of Chemical Engineering, Symposium Series, 142, 553–564. Boston, J. F., & Sullivan, S. L. (1974). A new class of solution methods for multicomponent, multistage separation processes. CanadianJournal of Chemical Engineering, 52, 52–63. Chang, Y. A., & Seader, J. D. (1988). Simulation of continuous reactive distillation by homotopy—Continuation method. Computers and Chemical Engineering, 12(12), 1243–1255. Chen, L., Hontoir, Y., Huang, D., Zhang, J., & Julian Morris, A. J. A. (2004). Combining first principles with black-box techniques for reaction systems. Control Engineering Practice, 12(7), 819–826. Doherty, M. F., & Buzad, G. (1992). Reactive distillation by design. Chemical Engineering Research and Design, Transactions of the Institute of Chemical Engineers, Part A, 70, 448–458.
239
Finlayson, B. A. (1980). Nonlinear analysis in chemical engineering. New York, NY: McGraw-Hill. Gill, P. E., Murray, W., & Saunders, M. A. (1988). User’s Guide for SNOPT 5.3/6.0: A Fortran Package for Large Scale Nonlinear Programming, Systems Optimization Laboratory, Stanford University, CA. Grosser, J. H., Doherty, M. F., & Malone, M. F. (1987). Modeling of reactive distillation systems. Industrial and Engineering Chemistry Research, 26, 983–989. Gru¨ner, S., Mohl, K.-D., Kienle, A., Gilles, E. D., Fernholz, G., & Friedrich, M. (2003). Nonlinear control of a reactive distillation column. Control Engineering Practice, 11(8), 915–925. Jacobs, R., & Krishna, R. (1993). Multiple solutions in reactive distillation for methyl tert-butyl ether synthesis. Industrial and Engineering Chemistry Research, 32, 1706–1709. Jones, D. I., & Finch, J. W. (1984). Comparison of optimization algorithms. International Journal of Control, 40, 747. Khaledi, R., & Young, B. R. (2005). Modeling and model predictive control of composition and conversion in an ETBE reactive distillation column. Industrial and Engineering Chemistry Research, 44, 3134–3145. Krishnamurthy, R., & Taylor, R. (1985). A nonequilibrium stage model of multicomponent separation processes. AIChE Journal, 32, 449–465. Kumar, A., & Daoutidis, P. (1999). Modeling, analysis and control of ethylene glycol reactive distillation column. AIChE Journal, 32, 449–465. Marlin, T. E. (1995). Process control: Design processes and control systems for dynamic performance. New York: McGraw Hill. Meadows, E. S., & Rawlings, J. B. (1997). Model predictive control. In M. Henson, & D. Seaborg (Eds.), Nonlinear process control. Prentice-Hall PTR. Riggs, J. B. (2001). Chemical process control (2nd ed.). Lubbock, TX: Ferret Publishing. Ruiz, C. A., Basualdo, M. S., & Scenna, N. J. (1995). Reactive distillation dynamic simulation. Chemical Engineering Research and Design, Transaction of the Institution of Chemical Engineering, Part A, 73, 363–378. Seferlis, P., & Grievink, J. (2001). Optimal design and sensitivity analysis of reactive distillation units using collocation models. Industrial and Engineering Chemistry Research, 40, 1673–1685. Sneesby, M. G., Tade, M. O., Datta, R., & Smith, T. N. (1997). ETBE synthesis via reactive distillation. 2. Dynamic simulation and control aspects. Industrial and Engineering Chemistry Research, 36, 1870–1881. Sneesby, M. G., Tade, M. O., & Smith, T. N. (1998). Multiplicity and pseudo-multiplicity in MTBE and ETBE reactive distillation. Chemical Engineering and Research Design, Transactions of the Institution of Chemical Engineers, Part A, 76, 525–531. Soave, G. (1992). Chemical and Engineering Science, 27, 1197. Suzuki, I., Yagi, H., Komatsu, H., & Hirata, M. (1971). Calculation of multicomponent distillation accompanied by chemical reaction. Journal of Chemical Engineering of Japan, 4, 26–33. Tyreus, B. D., & Luyben, W. L. (1992). Tuning PI controllers for integrator/deadtime processes. I&ECR, 31, 2625–2635. Vo¨lker, M., Sonntag, C., & Engell, S. (2006). Control of integrated processes: A case study on reactive distillation in a medium-scale pilot plant, Control Engineering Practice, in press, Corrected Proof, Available online 5 June 2006, doi:10.1016/j.conengprac.2006.03.00.
Control Engineering Practice 15 (2007) 231–239 www.elsevier.com/locate/conengprac
Nonlinear model predictive control of a reactive distillation column Rohit Kawathekar, James B. Riggs Chemical Engineering, Texas Tech University, Lubbock, TX 79409-3121, USA Received 5 July 2005; accepted 10 July 2006 Available online 7 September 2006
Abstract This paper considers the application of nonlinear model predictive control (NLMPC) to a highly nonlinear reactive distillation column. NLMPC was applied as a nonlinear programming problem using orthogonal collocation on finite elements to approximate the ODEs that constitute the model equations for the reactive distillation column. Diagonal PI controls were used to identify that the [L/D,V] and the [L/D,V/B] configurations performed best. NLMPC was applied using the [L/D,V] configuration and found to provide a factor of 2–3 better performance than the corresponding PI controller. The effect of process/model mismatch on the performance of the NLMPC controller was also evaluated. r 2006 Elsevier Ltd. All rights reserved. Keywords: Reactive distillation; Nonlinear model predictive control; Configuration selection; Process/model mismatch
1. Introduction Reactive distillation combines both separation and reaction in one unit and has been applied industrially for a number of years. Reactive distillation can offer significant economic advantages for certain cases, particularly for systems that involve reversible reactions. Examples of commercial successes of reactive distillation include nylon 6,6 processes, methyl acetate processes, ethyl acetate processes and methyl tert-butyl ether processes (Doherty & Buzad, 1992). Most of the literature available on reactive distillation is based on steady-state conditions including process design (e.g., Chang & Seader, 1988; Krishnamurthy & Taylor, 1985) and the analysis of multiple steady states (e.g., Abufares & Douglas, 1995; Jacobs & Krishna, 1993; Sneesby, Tade, & Smith, 1998). Dynamic modeling and simulation (e.g., Grosser, Doherty, & Malone, 1987; Ruiz, Basualdo, & Scenna, 1995) have also been studied, but a relatively small amount of work has been reported on the control of reactive distillation columns. Chen, Hontoir, Huang, Zhang, and Julian Morris (2004) performed a nonlinear dynamic modeling case study of a reactive Corresponding author. Tel.: +1 806 742 1765; fax: +1 806 742 3552.
E-mail address: [email protected] (J.B. Riggs). 0967-0661/$ - see front matter r 2006 Elsevier Ltd. All rights reserved. doi:10.1016/j.conengprac.2006.07.004
distillation column using a nonlinear first-principles model combined with an artificial neural network model. Reactive distillation is a challenge for control due to process nonlinearity and complex interactions between vapor–liquid equilibrium and chemical reactions. Al-Arfaj and Luyben (2000), Sneesby, Tade, Datta, and Smith (1997) and Kumar and Daoutidis (1999) have considered decentralized PI control structures for reactive distillation columns. Khaledi and Young (2005) applied PI controls and linear model predictive control to a dynamic simulation of an ETBE reactive distillation column using the (L,V) configuration. The linear model predictive controller used first-order plus deadtime models for each of its input/output models. Kumar and Daoutidis (1999) have discussed the superior performance of nonlinear controllers compared to linear controllers for reactive distillation systems. Vo¨lker, Sonntag, and Engell (2006) applied a linear reduced-order MPC controller to a semibatch reactive distillation column. Gru¨ner et al. (2003) developed a general controller for reactive distillation systems based on asymptotically exact input/output linearization and applied their controller to a simulation of an industrial reactive distillation column with improved performance compared to a well-tuned linear controller. This paper considers the application of nonlinear model predictive control (NLMPC) to an ethyl acetate reactive
ARTICLE IN PRESS R. Kawathekar, J.B. Riggs / Control Engineering Practice 15 (2007) 231–239
232
Nomenclature
xi,j
B D Fi Hi hi ki
yi,j
kfi L L/D Mi l nc Pi Psi,j QC QR QMS ri u V V/B Voli
bottom product flow rate distillate product flow rate the external feed rate to the ith stage. the enthalpy of the vapor leaving the ith stage the enthalpy of the liquid leaving the ith stage the reaction equilibrium rate constant for the ith stage the reaction rate constant for the ith stage the reflux rate the reflux ratio the molar holdup of liquid on the ith stage the number of components the pressure on the ith stage the vapor pressure of component j on the ith stage the condenser duty the reboiler duty the mover suppression factor for the NLMPC controller the reaction rate on the ith stage the manipulate variable the boilup rate the boilup ratio the volume of the liquid holdup on the ith stage
distillation column. In Section 2, the case study is described in detail. Application of diagonal PI controls is presented in Section 3. Development and implementation of NLMPC algorithm is discussed in Section 4. In Section 5, the comparison between PI and NLMPC control results is presented. The effect of process-model mismatch on the closed-loop performance of NLMPC controller is analyzed in Section 6. 2. Model development This study is based on the production of ethyl acetate by a reactive distillation using the esterification reaction between acetic acid and ethanol. The achievable conversion in this reversible reaction is limited by the equilibrium conversion. The reaction is slightly endothermic and takes place in the liquid phase. Though the esterification reaction is self-catalyzed, sulfuric acid can act as external catalyst to enhance the reaction rate: CH3 COOH þ C2 H5 OH Sulfuric acid
!
CH3 COOC2 H5 þ H2 O:
ð1Þ
The reaction kinetics for this reaction were proposed by Aljeski and Duprat (1996) and include reaction equilibrium limitations. Bock and Wozny (1997) reported a detailed analysis of ethyl acetate reactive distillation column. The slow reaction rate for this system leads to lower conversion
ysp
the mole fraction of the jth component in the liquid on the ith stage the mole fraction of the jth component in the vapor on the ith stage the setpoint for the controlled variable
Greek Symbols gi,j mj ri
the liquid phase activity coefficient for the jth component on the ith stage the stochiometric coefficient for the jth component in the primary reaction the liquid density for the ith stage
Subscripts A B C c D i f j n
acetic acid ethyl alcohol ethyl acetate condenser water stage number feed component number the top tray number
for the reactive distillation column. A reactive column with countercurrent flow of reactants was used. Pure ethanol and acetic acid were fed separately into a column that operates at atmospheric pressure. Under these conditions, acetic acid is the heaviest of the components and moves toward the bottom of the column. Ethyl acetate is the lightest and moves toward the top of the column. It is expected that the middle portion of the distillation column is the primary reaction zone. The rectifying section fractionates the ethyl acetate out of acetic acid, and the stripping section removes alcohol from water. This quaternary system consisting of ethanol, acetic acid, water and ethyl acetate is highly nonideal. It can form four binary azeotrope mixtures and one ternary azeotrope. Over a wide range of composition, ethanol and water do not differ greatly in volatility, making it difficult to produce only water as the bottom product. The reactant ethanol has a relatively high volatility in the reaction zone. This leads to a low composition of ethanol in the liquid phase, reducing the production rate of ethyl acetate. Under these conditions, the rate of the esterification reaction between acetic acid and ethanol is generally low, which implies that it is favored by long residence times in each stage. It is evident from all the previous studies on ethyl acetate reactive distillation columns that an unfavorable physical equilibrium makes the production of highpurity ethyl acetate impossible from a single distillation column.
ARTICLE IN PRESS R. Kawathekar, J.B. Riggs / Control Engineering Practice 15 (2007) 231–239 Reactive Column
AA feed
Recovery Column
101.3 kPa
350.3 kPa Purge
EtOH feed
RX
RC
Ethyl acetate product Recycle
Fig. 1. Two-column sequence for the production of high-purity ethyl acetate.
However, Seferlis and Grievink (2001) proposed adding a recovery column (Fig. 1), which operates at a higher pressure to produce a high-purity ethyl acetate product. The recovery column was combined with the reactive distillation column and designed to operate at a higher pressure (350 kPa) to break the azeotropes, which increases the overall conversion and produces a high-purity ethyl acetate product. The distillate stream from the reactive column is fed to the recovery column. At a relatively high pressure, ethyl acetate becomes heavier than ethanol and water so that it appears as the bottom product. The target purity level of the ethyl acetate was set at 99.5%. Reaction in the recovery column is negligible because the column operates without any sulfuric acid catalyst and the liquid holdup on the trays is small. Virtually all of the acetic acid in the overhead of the reactive column ends up in the bottom product of the recovery column, and therefore, directly affects the purity of final ethyl acetate product. Hence, a specification is imposed on the maximum allowable concentration of acetic acid in the distillate of the reactive column. To assist in the reduction of the acetic acid in the distillate of the reactive column, the reactive column was designed with a reduced liquid phase holdup in the upper section of trays of the reactive column to suppress the reverse reaction of ethyl acetate. Seferlis and Grievink (2001) developed an economically optimized steady-state design for a reactive distillation process based on the two-column configuration (Fig. 1). The exact design and steady-state operating conditions of the reactive distillation column from Seferlis and Grievink (2001) were used as the basis of this study while the recovery column was not directly modeled. An equilibrium stage dynamic model of the reactive distillation column was developed and benchmarked against the results
233
reported by Seferlis and Grievink (2001). This dynamic model served as the process that was controlled by the PI and NLMPC controllers in this study. An equilibrium stage dynamic model of the reactive distillation column was implemented with the following assumptions: (1) each tray was assumed to be an ideal stage, (2) nonequal molar overflow was considered, (3) the liquid activity coefficient was modeled using the approach of Suzuki, Yagi, Komatsu, and Hirata (1971) [15], (4) the SRK method (Soave, 1992) was used to model the enthalpy departure functions, (5) a hydraulic time constant was used to model the liquid flow from each tray, (6) tray temperatures in the rectifying and stripping sections were used to infer the impurity composition of the overhead and bottom products, respectively, (7) the heat of dilution of acetic acid in water was neglected and (8) equilibrium limited reaction kinetics (Aljeski & Duprat, 1996) are assumed. For simulation purposes, the process is modeled using a detailed tray-to-tray model. For a standard reactive distillation tray where vapor holdup is considered negligible and reaction takes place in liquid phase, the modeling equations can be written as follows: nc X dðM li Þ ¼ F i þ Liþ1 þ V i1 Li V i þ Vol i m j ri , dt j¼1
dðM li xi;j Þ ¼ F i xi;f þ Liþ1 xiþ1;j þ V i1 yi1;j dt Li xi;j V i yi;j þ Vol i mj ri ,
(2)
ð3Þ
dðM li hi Þ ¼ F i hf i þ Liþ1 hiþ1 V i1 H i1 Li hi V i H i , dt
(4)
dhi 0, dt
(5)
1 ri ¼ r2i kfi xA;i xB;i xC;i xD;i , ki
(6)
yi;j Pi ¼ gi;j xi;j Psi;j .
(7)
For the partial condenser, the material balance and energy balance equations are as follows: nc X dðM lc Þ ¼ V n Lc V D þ Vol c m j rc , dt j¼1
(8)
dðM lc xc;j Þ ¼ V n yn;j Lc xc;j V D yD;j þ Vol c mj rc , dt
(9)
dðM lc hc Þ ¼ V n H n Lc hc V D H D Qc , dt
(10)
dhlc 0, dt
(11)
ARTICLE IN PRESS R. Kawathekar, J.B. Riggs / Control Engineering Practice 15 (2007) 231–239
234
ri ¼
r2i kfi
1 xA;i xB;i xC;i xD;i , ki
(12)
yi;j Pi ¼ gi;j xi;j Psi;j ,
(13)
where subscript c represents condenser stage, subscript n represents the top tray of the column and subscript D represents the distillate. VD represents the vapor distillate from the partial condenser stage. The reflux accumulator stage is modeled as a liquid holdup tank. The terms relating to the vapor phase, i.e., vapor flow, vapor composition, vapor enthalpy corresponding to reflux accumulator are absent in the material and energy balances. The material and energy balance equations for the reboiler are written as follows: nc X dðM lB Þ ¼ Li B V B þ Vol B mj rB , dt j¼1
(14)
dðM lB xB;j Þ ¼ Li xi;j BxB;j V B yB;j þ Vol B mj rB , dt
(15)
dðM lB hB Þ ¼ Li hi BhB V B H B þ QR , dt
(16)
dhlB 0, dt
(17)
1 ri ¼ r2i kfi xA;i xB;i xC;i xD;i , ki
(18)
yi;j Pi ¼ gi;j xi;j Psi;j .
(19)
D
PC
LT
RSP LC
CW LC
LT PT
L
Fig. 2. Schematic of the column pressure controls and the accumulator level controls.
The overhead level and pressure control structures are shown in Fig. 2 based on producing a vapor distillate product. The overhead pressure control loop sets the distillate vapor flow rate. The overhead pressure was modeled considering the total number of moles of vapor in the overhead and the vapor volume in the overhead line and the overhead condenser. A flooded condenser is used and the level of liquid in the condenser is used to vary the surface area in the condenser available for condensing the overhead vapor. An accumulator level controller determines the setpoint for the condenser level controller. Because all the acetic acid that leaves in the overhead of the reactive distillation column ends up in the ethyl acetate product produced by the recovery column, the primary composition control point is the acetic acid concentration in the overhead product from the reactive distillation column (i.e., the design setpoint is 0.04% acetic acid in the overhead). The composition setpoint for the bottoms is 0.08% ethyl acetate, which was determined by Seferlis and Grievink (2001) using a steady-state economic analysis of the two-column system. The temperature of the second tray from the top was used to infer the impurity level in the overhead product and the fifth tray from the bottom was used to infer the impurity level in the bottom product. These tray temperatures were identified using a steady-state model of the reactive distillation column to determine the tray temperatures that exhibited the strongest correlation between changes in the product purities and changes in the tray temperatures. The correlations between tray temperature and product composition were updated using new composition analyzer information when it became available. It was determined that this inferential tray temperature scheme provided excellent estimates of the product composition by comparing the inferential estimate of the product impurity levels with the corresponding values generated by the process simulator. Therefore, an observer was not considered necessary in this case. Table 1 lists the steady-state design conditions for the reactive distillation column studied here. The model equations for the reactive distillation column (i.e., the liquid dynamics equation, the component material balances and the energy balance for each tray, the component material balances and the overall material balances for the accumulator, reboiler and condenser, and the PI control equations for the pressure controller and the level controllers) were integrated numerically using the Euler method with a step size of 0.5 s. To improve the computational efficiency of the dynamic reactive distillation model, the inside-out algorithm (Boston & Sullivan, 1974) was used to calculate the component K-values. If a tray temperature changed by more than 1 1C, the parameters of the K-value correlation were updated. In addition, all K-value correlations were updated each 30 s of simulated time. The accumulator and reboiler level controllers were tuned for critically damped behavior.
ARTICLE IN PRESS R. Kawathekar, J.B. Riggs / Control Engineering Practice 15 (2007) 231–239
4.6E-04
31 25–3a 7–24 1–6
Stage holdups Rectifying section (m3) Reactive section (m3) Stripping section (m3)
0.151 1.75 1.65
Distillate flow (kmol/s) Bottoms flow (kmol/s) Reflux ratio
4.88 103 9.09 104 0.8211
Feed flow rate Acetic acid feed (kmol/s) Ethanol feed (kmol/s) Recycle feed (kmol/s)
8.92 104 8.74 104 4.03 103
Distillate composition (mol frac) Ethanol Acetic acid Water Ethyl acetate
0.2383 0.0004 0.2093 0.5520 0.0249 0.0876 0.8867 0.0008 171.0
Condenser duty (kW)
137.7
Low holdup was used to reduce reaction in the rectifying section.
Table 2 Steady-state changes for a 10% change in a manipulated variable Configuration
[L/D,V] [L/D,V/B] [L,V] [L,V/B]
L/D,V L/D,V/B L/D,B
4.4E-04 4.2E-04 4.0E-04 3.8E-04 0
500
1000
1500
2000
Time (min) Fig. 3. Overhead composition control using PI controls for an unmeasured feed composition upset.
3. Diagonal PI control results
Bottoms composition (mol frac) Ethanol Acetic acid Water Ethyl acetate Reboiler duty (kW)
a
Overhead impurity (mole fraction)
Table 1 Steady-state design conditions for the reactive distillation column Number of stages Rectifying section Reactive section Stripping section
235
% change in the process gain Overhead
Bottom
200 150 150 33
300 300 600 300
The steady-state gains for the ethyl acetate reactive distillation column are listed in Table 2. Note that for approximately a 10% change in each manipulated variable of each configuration, the process gain mostly changed by 150–600%, indicating severe steady-state process nonlinearity. In addition, several of the off-diagonal gains demonstrated bi-directionality. As is usually the case, process integration (i.e., adding the recovery column with recycle in this case) offers significant steady-state economic benefits but at the expense of a much more challenging control problem. The dynamic model of the reactive distillation column described in this section was used as the process, which was controlled by PI and nonlinear MPC controllers.
Diagonal PI composition controllers were used to compare the control performance of various control configurations and as a benchmark for the NLMPC control results. Because the primary control objective for the reactive distillation column is control of the impurity in the overhead product, the bottom composition control loop was tuned first for a closed-loop overdamped response. Then the overhead PI composition controller was tuned using setpoint changes. The overhead tuning procedure was based on using an ATV test using a single tuning factor (Riggs, 2001) combined with Tyreus– Luyben settings (Tyreus & Luyben, 1992) and finally adjusted (Riggs, 2001) to minimize the IAE from setpoint for a series of setpoint changes. The control performance of the overhead and bottoms composition controllers was evaluated for several unmeasured disturbances in feed rate and feed composition. Because fresh acetic acid and ethyl alcohol are each fed to the reactive column as well as the recycle stream from the recovery column, feed composition upsets were implemented by changing the ratio of fresh acetic acid and ethyl alcohol and feed rate changes were implemented by changing the recycle feed rate. The [L,B], [L,V], [L,V/B], [L/D,B], [L/D,V] and [L/D,V/ B] configurations were compared for control performance based on dual composition PI control (note that for a [A,B] configuration A is manipulated to control the overhead composition and B to control the bottom composition). Configurations using the distillate rate, D, as a manipulated variable were not considered because D is used to control the column pressure (Fig. 2). Figs. 3 and 4 show the overhead and bottoms composition control performance for the [L/D,B], [L/D,V] and [L/ D,V/B] configurations for an unmeasured feed composition upset. Table 3 lists the IAEs for each configuration for a feed rate upset. Based on these results, the [L/D,V] and [L/ D,V/B] configurations provided the best overall control performance.
ARTICLE IN PRESS R. Kawathekar, J.B. Riggs / Control Engineering Practice 15 (2007) 231–239
236
differences observed between the process output and the model prediction are due to additive step disturbances in the output. These disturbance terms are also assumed to remain constant over the prediction horizon.
Bottom ethyl acetate impurity (mole fraction)
8.3E-04 8.2E-04 8.1E-04 8.0E-04 7.9E-04
4.1. Solution algorithm
7.8E-04 7.7E-04
L/D,V/B
7.6E-04
L/D,V
7.5E-04
L/D,B
7.4E-04 0
2000
4000
6000
Time (min)
Fig. 4. Bottom composition control performance using PI control for an unmeasured feed composition upset.
Table 3 Dual PI composition control performance indices for unmeasured feed rate disturbance Configuration
IAE
[L,V] [L,B] [L,V/B] [L/D,V] [L/D,B] [L/D,V/B]
Overhead
Bottoms
0.143 0.177 0.142 0.115 0.393 0.114
1.278 0.830 1.239 1.270 1.215 0.988
4. Implementation of NLMPC NLMPC uses nonlinear models in a model predictive control framework to choose control action. The NLMPC controller chooses the future control moves that minimize the following objective function: min
N 1 X j¼0
½ysp yj 2 þ
N 1 X
QMS ½Duj 2 ,
(20)
j¼0
subject to bounds on manipulated, output and state variables. The first term represents the error from setpoint and the second term represents the weighted variation in the manipulated variables. Once the optimum set of future moves is determined, the initial control move is implemented. In the next control cycle, the next set of optimum control moves is determined. Even though the full set of optimum control moves into the future are calculated at each control interval, only the first move is actually implemented at each control interval. The formation of the NLMPC objective function will involve a disturbance term, which is added to the output prediction over the entire prediction horizon to match the model prediction to the current process measurement. The NLMPC procedure in this work assumes that the
NLMPC was implemented by setting up the control problem as a nonlinear programming (NLP) problem and solving it over the prediction horizon. It is necessary to simultaneously solve an optimization problem (the tradeoff between control to setpoint and changes in the manipulated variables) and the system model equations. These two procedures may be implemented either sequentially or simultaneously. 4.1.1. Sequential solution and optimization algorithm The sequential algorithm employs separate algorithms to solve the differential equations and to carry out the optimization. First, a manipulated variable profile is selected and the differential equations are solved numerically to obtain the controlled variable profile. The objective function is then determined. The gradient of the objective function with respect to the manipulated variables is determined either by numerical perturbation or by using analytical derivatives. The control profile is then updated using an optimization algorithm. The process is repeated until the optimal profile is obtained. This is referred as a sequential solution and optimization algorithm. The implementation of the sequential solution and optimization algorithm is relatively simple to apply. However, there are some drawbacks associated with this approach. The sequential solution and optimization requires the solution of differential equations at each iteration of the optimization. Jones and Finch (1984) found that such methods spend about 85% of the time integrating the model equations in order to obtain gradient information. This can make the implementation of this algorithm computationally expensive for cases involving a large number of model equations. The gradient information required for the optimization procedure is often obtained through numerical differentiation, as the analytical derivatives are not always available for highly nonlinear model equations involving complicated thermodynamic relations. To obtain the gradients using finite difference typically involves differencing the output of an integration routine with adaptive step sizes. Gill, Murray, and Saunders (1988) points out that the integration error is unpredictable and hence differencing output of an integration routine greatly degrades the quality of the finite difference derivatives. It is also difficult to incorporate the constraints on state variables with the use of the sequential solution and optimization approach (Meadows & Rawlings, 1997). 4.1.2. Simultaneous solution and optimization algorithm The simultaneous solution and optimization algorithm involves the model equations appended to the optimization
ARTICLE IN PRESS R. Kawathekar, J.B. Riggs / Control Engineering Practice 15 (2007) 231–239
problem as equality constraints. Then the NLP problem is posed to optimize the objective function such that the algebraic form of the model differential equations are satisfied, and the upper and lower constraints on states, controlled variables and manipulated variables are met. This can greatly increase the size of optimization problem, leading to a tradeoff between the two approaches. Meadows and Rawlings (1997) reported that for small problems with few states and a short prediction horizon, the sequential solution and optimization algorithm is more computationally efficient. For larger problems, the simultaneous solution and optimization approach is more computationally efficient and more reliable. In the simultaneous solution and optimization approach, the model equations and the process constraints (e.g., a manipulated variable cannot be less than zero) become constraints for the optimization problem to determine the future control moves. Since the model equations largely appear as ODEs, orthogonal collocation on finite elements (Finlayson, 1980) was used to convert the ODEs in time into a set of algebraic equality constraints for the optimization problem. Due to the size of the resulting NLMPC problem considered here, the simultaneous solution approach was chosen for this work. The process model and the NLMPC controller were implemented in Fortran code. The NLP problem resulting from the implementation of NLMPC was solved using the SNOPT software (Gill et al., 1988). Typical computational times for 30 h of simulation of the NLMPC controller using an 1100 MHz in an AMD PC with Windows XP required 10 h CPU. 4.2. Selection of tuning parameters The tuning parameters that have a significant effect on NLMPC performance are the prediction horizon, control horizon, sampling interval, penalty weight matrices and move suppression factors. A set of heuristics based on the linear systems and numerical simulations were used to select the final tuning parameters Following is a discussion of the tuning parameters used for the tuning of each NLMPC controller considered: (1) The control interval was set at 20 min, which was shown to be equivalent to continuous control using a correlation by Marlin (1995). (2) The prediction horizon was set at 80 by using different levels and comparing control performance. (3) A control horizon of 15 was found to provide good control performance. (4) The equal concern errors were selected to provide a 10 times greater weighting of the overhead errors from setpoint than the bottoms because impurities in the overhead directly affect the purity of the ethyl acetate final product. (5) The move suppression factors were adjusted to minimize the IAE from setpoint for the overhead product for a series of overhead impurity setpoint changes. 5. Comparison between PI and NLMPC The NLMPC control algorithm described in the previous sections was applied for dual-ended composition control of
237
Table 4 Ethyl acetate NLMPC control performance indices for overhead impurity setpoint tracking Configuration
Overhead loop IAE
Bottoms loop IAE
[L/D,V] NLMPC [L/D,V] PI
0.91 2.41
2.14 6.84
Table 5 NLMPC control performance indices for unmeasured feed composition disturbance rejection Configuration
[L/D,V] NLMPC [L/D,V] PI
IAE Overhead
Bottoms
0.0609 0.115
0.483 1.270
the simulation of the ethyl acetate reactive distillation column. The previous PI dual-ended composition control results indicated that [L/D, V] configuration provided good control performance for setpoint tracking as well as the rejection of unmeasured disturbances. NLMPC was applied using the [L/D, V] configuration to determine the benefits compared to conventional PI controls. Table 4 summarizes the comparisons between diagonal PI control and NLMPC using the [L/D,V] configuration for a series of setpoint changes in the overhead impurity. Table 5 summarizes the comparisons between diagonal PI control and NLMPC for an unmeasured feed composition upset. These results are based on an NLMPC controller with an almost perfect process model. Based on the IAEs of these tests, the NLMPC controller outperformed the PI controller by a factor of 2–3. Due to the severe nonlinearity of this process, one would expect a larger disparity between PI and NLMPC. On the other hand, this column was operated largely as a single-ended column control problem because the bottom composition control was detuned compared to the more important overhead product, which significantly simplified the control problem. 6. Process/model mismatch Process/model mismatch for the NLMPC controller was represented by introducing a 5% and a 25% error between the reaction equilibrium constant used in the dynamic process simulator and the model equations used by the NLMPC controller to calculate control action. Figs. 5–9 show the effect of process model/mismatch for NLMPC for an unmeasured feed composition upset. The PI control results are also shown for this case. To more completely analyze the sensitivity of PI and NLMPC control performance to unmeasured disturbances, a sinusoidal disturbance in the recycle feed rate with constant amplitude and varying frequency was applied to the process with a PI controller and an NLMPC controller
ARTICLE IN PRESS R. Kawathekar, J.B. Riggs / Control Engineering Practice 15 (2007) 231–239
238
1.68E+02 Reboiler duty (btu/sec)
Overhead impurity (mole fraction)
4.2E-04 4.0E-04 3.8E-04 0% mismatch 5% mismatch 25% mismatch PI
3.6E-04 3.4E-04 0
200
400
600
Time (min)
25% mismatch PI
1.62E+02 1.60E+02 1.58E+02 1.56E+02 1.54E+02 500
1500
0.07 PI
Overhead impurity amplitude ratio
0.06
8.1E-04
7.9E-04 0% mismatch 5% mismatch 25% mismatch PI
7.7E-04
NLMPC
0.05
NLMPC
0.04
5% mismatch
0.03 0.02 0.01 0 0
7.5E-04 0
500
1000
Time (min)
Fig. 6. Control results for the bottom product for an unmeasured feed composition upset showing the effect of process/model mismatch for NLMPC and compared with PI control for the [L/D,V] configuration.
1.3 1.2 1.1 1.0 0% mismatch 0.9
5% mismatch 25% mismatch
0.8
PI 0.7 0
1000 Time (min)
Fig. 8. The manipulated variable profile for bottoms for an unmeasured feed composition upset for the [L/D,V] configuration.
8.3E-04 Bottoms ethyl acetate impurity (mole fraction)
5% mismatch
1.64E+02
0
Fig. 5. Control results for the overhead product for an unmeasured feed composition upset showing the effect of process/model mismatch for NLMPC and compared with PI control for the [L/D,V] configuration.
Reflux ratio
0% mismatch
1.66E+02
200
400
600
Time (min)
Fig. 7. The manipulated variable profile for overhead for an unmeasured feed composition upset for the [L/D,V] configuration.
in use. A sinusoidal disturbance to the recycle feed rate was applied with an amplitude of 10% of the base case recycle feed rate. The closed-loop amplitude ratio for the overhead impurity was calculated by taking the ratio of the amplitude of overhead impurity resulting from the
500
1000
1500
Time period for sinusoidal feed disturbance (min)
Fig. 9. A comparison of closed-loop amplitude ratios for PI and NLMPC for the [L/D,V] configuration.
sinusoidal disturbance to the amplitude of the recycle feed disturbance. The closed-loop performances of the PI control structure, NLMPC controller as well as NLMPC with 5% process/model mismatch were compared using the unmeasured sinusoidal recycle feed rate disturbance. A plot of the closed-loop amplitude ratio versus the period of the sinusoidal feed disturbance (Fig. 9) shows that for the PI controller as well as the NLMPC controller exhibit sensitivity to the feed disturbance in a specific range of frequencies (120 min period to a 720 min period). The closed-loop amplitude ratio for the diagonal PI control structure is shown to be 2–3 times larger than that for NLMPC over the sensitive range. At low periods (high frequencies), the inertia of the process damps out the effect of the sinusoidal disturbance. At large periods (low frequencies), the feedback controllers have ample time to absorb the disturbances by feedback. Fig. 9 also indicates that the 5% process/model mismatch does not significantly affect the control performance for the NLMPC controller. These results are consistent with the previous results for the comparisons between PI and NLMPC control performance for the setpoint tracking as well as unmeasured disturbance rejection, but more quantitatively define the
ARTICLE IN PRESS R. Kawathekar, J.B. Riggs / Control Engineering Practice 15 (2007) 231–239
control performance differences between PI and NLMPC for this case. 7. Conclusions As expected, NLMPC was found to provide significantly better control performance than PI controls, but this work quantifies the advantage as a factor of 2–3 reduction in variability. The advantage of the NLMPC controller comes from faster closed-loop dynamic performance compared to the PI controller resulting from using a nonlinear dynamic model of the reactive distillation column considered here. In addition, NLMPC was found not to be particularly sensitive to process/model mismatch. Even though the reactive distillation column considered in this study was shown to be extremely nonlinear, the PI control was able to control the process reasonably well. Acknowledgements This work was supported by the member companies of the Texas Tech Process Control and Optimization Consortium. The authors would also like to acknowledge assistance from Charles R. Cutler, Cutler-Technologies. References Abufares, A. A., & Douglas, P. L. (1995). Mathematical modeling and simulation of an MTBE catalytic distillation process using SPEEDUP and aspen plus. Chemical Engineering Research and Design, Transactions of the Institute of Chemical Engineers Part A, 73, 3–12. Al-Arfaj, M. A., & Luyben, W. L. (2000). Comparative control of ideal and methyl acetate reactive distillation. Chemical Engineering Science, 57, 5039–5050. Aljeski, K., & Duprat, F. (1996). Dynamic simulation of the multicomponent reactive distillation. Chemical Engineering Science, 51, 4237–4252. Bock, H., & Wozny, G. (1997). Analysis of distillation and reaction rate in reactive distillation. Distillation and absorption ‘97. Institute of Chemical Engineering, Symposium Series, 142, 553–564. Boston, J. F., & Sullivan, S. L. (1974). A new class of solution methods for multicomponent, multistage separation processes. CanadianJournal of Chemical Engineering, 52, 52–63. Chang, Y. A., & Seader, J. D. (1988). Simulation of continuous reactive distillation by homotopy—Continuation method. Computers and Chemical Engineering, 12(12), 1243–1255. Chen, L., Hontoir, Y., Huang, D., Zhang, J., & Julian Morris, A. J. A. (2004). Combining first principles with black-box techniques for reaction systems. Control Engineering Practice, 12(7), 819–826. Doherty, M. F., & Buzad, G. (1992). Reactive distillation by design. Chemical Engineering Research and Design, Transactions of the Institute of Chemical Engineers, Part A, 70, 448–458.
239
Finlayson, B. A. (1980). Nonlinear analysis in chemical engineering. New York, NY: McGraw-Hill. Gill, P. E., Murray, W., & Saunders, M. A. (1988). User’s Guide for SNOPT 5.3/6.0: A Fortran Package for Large Scale Nonlinear Programming, Systems Optimization Laboratory, Stanford University, CA. Grosser, J. H., Doherty, M. F., & Malone, M. F. (1987). Modeling of reactive distillation systems. Industrial and Engineering Chemistry Research, 26, 983–989. Gru¨ner, S., Mohl, K.-D., Kienle, A., Gilles, E. D., Fernholz, G., & Friedrich, M. (2003). Nonlinear control of a reactive distillation column. Control Engineering Practice, 11(8), 915–925. Jacobs, R., & Krishna, R. (1993). Multiple solutions in reactive distillation for methyl tert-butyl ether synthesis. Industrial and Engineering Chemistry Research, 32, 1706–1709. Jones, D. I., & Finch, J. W. (1984). Comparison of optimization algorithms. International Journal of Control, 40, 747. Khaledi, R., & Young, B. R. (2005). Modeling and model predictive control of composition and conversion in an ETBE reactive distillation column. Industrial and Engineering Chemistry Research, 44, 3134–3145. Krishnamurthy, R., & Taylor, R. (1985). A nonequilibrium stage model of multicomponent separation processes. AIChE Journal, 32, 449–465. Kumar, A., & Daoutidis, P. (1999). Modeling, analysis and control of ethylene glycol reactive distillation column. AIChE Journal, 32, 449–465. Marlin, T. E. (1995). Process control: Design processes and control systems for dynamic performance. New York: McGraw Hill. Meadows, E. S., & Rawlings, J. B. (1997). Model predictive control. In M. Henson, & D. Seaborg (Eds.), Nonlinear process control. Prentice-Hall PTR. Riggs, J. B. (2001). Chemical process control (2nd ed.). Lubbock, TX: Ferret Publishing. Ruiz, C. A., Basualdo, M. S., & Scenna, N. J. (1995). Reactive distillation dynamic simulation. Chemical Engineering Research and Design, Transaction of the Institution of Chemical Engineering, Part A, 73, 363–378. Seferlis, P., & Grievink, J. (2001). Optimal design and sensitivity analysis of reactive distillation units using collocation models. Industrial and Engineering Chemistry Research, 40, 1673–1685. Sneesby, M. G., Tade, M. O., Datta, R., & Smith, T. N. (1997). ETBE synthesis via reactive distillation. 2. Dynamic simulation and control aspects. Industrial and Engineering Chemistry Research, 36, 1870–1881. Sneesby, M. G., Tade, M. O., & Smith, T. N. (1998). Multiplicity and pseudo-multiplicity in MTBE and ETBE reactive distillation. Chemical Engineering and Research Design, Transactions of the Institution of Chemical Engineers, Part A, 76, 525–531. Soave, G. (1992). Chemical and Engineering Science, 27, 1197. Suzuki, I., Yagi, H., Komatsu, H., & Hirata, M. (1971). Calculation of multicomponent distillation accompanied by chemical reaction. Journal of Chemical Engineering of Japan, 4, 26–33. Tyreus, B. D., & Luyben, W. L. (1992). Tuning PI controllers for integrator/deadtime processes. I&ECR, 31, 2625–2635. Vo¨lker, M., Sonntag, C., & Engell, S. (2006). Control of integrated processes: A case study on reactive distillation in a medium-scale pilot plant, Control Engineering Practice, in press, Corrected Proof, Available online 5 June 2006, doi:10.1016/j.conengprac.2006.03.00.