3528
Ind. Eng. Chem. Res. 2004, 43, 3528-3538
Dynamic Precompensation and Output Feedback Control of Integrated Process Networks Marie-Nathalie Contou-Carrere, Michael Baldea, and Prodromos Daoutidis* Department of Chemical Engineering and Materials Science, University of Minnesota, 421 Washington Avenue SE, Minneapolis, Minnesota 55455
In integrated process networks, the presence of large recycle flow rates induces a time scale separation where the individual units evolve in a fast time scale while the overall network evolves in a slow time scale. The slow dynamics of such networks is modeled by a high-index differentialalgebraic equation (DAE) system which, in the case of cascaded control configurations, has a control-dependent state space. In this paper, a minimal-order dynamic extension is proposed to obtain a modified DAE system of index 2 with a control-independent state space that can be subsequently used as the basis for output feedback controller design. The application of this method is illustrated for a reactor-condenser network and for a two-point control problem in a high-purity distillation column. 1. Introduction Process networks, consisting of reaction and separation units interconnected through material and energy recycle, constitute the rule rather than the exception in the process industries. The efficient transient operation of such networks is becoming increasingly important because the current environment of frequent changes in market conditions dictates frequent changes in operating conditions and targets.1,2 At the same time, recycle streams induce strong interactions between the individual process units, which can lead to complex overall “network” dynamics, including highly nonlinear behavior,3,4 inverse response,5 and even open-loop instability.6 Despite this fact, control methods that directly account for the nonlinear characteristics induced by recycle are really sparse. Most of the existing literature on the control of integrated process networks focuses on the design of linear multiloop control systems (e.g., refs 7 and 8) or addresses the problems of control structure selection and controller tuning within a broader plant-wide control framework (e.g., refs 9-12). In a different vein, a formal framework for stability analysis and stabilization of networks based on passivity and concepts from thermodynamics has recently been postulated.13,14 In this work, we address the nonlinear output feedback control of process networks with large material recycle compared to throughput. Our previous work15 has established that, owing to the presence of large and small flow rates, the dynamics of such networks exhibit a time scale separation, where the individual units evolve on a fast time scale and the overall dynamics on a slow time scale. The slow dynamics were shown to be modeled by a high-index differential-algebraic equation (DAE) system. To account for the time scale separation, we proposed a control strategy that involves distributed control for individual units in the fast time scale and supervisory control for the overall network in the slow time scale. In this framework, the design of nonlinear supervisory controllers on the basis of the DAE model * To whom correspondence should be addressed. Fax: (612) 626-7246. E-mail:
[email protected].
becomes important. This control problem is addressed through the derivation of an equivalent ordinary differential equation (ODE) representation of the slow dynamics, on the basis of which a controller is designed by application of well-known ODE control methods. Very often, cascaded control configurations are used in order to achieve a tighter coordination between distributed and supervisory controllers. Such an approach becomes necessary when the number of control objectives on the slow time scale exceeds the number of available manipulated inputs. In such cases, the constrained state space of the DAE model of the slow dynamics depends explicitly on the manipulated inputs (the slowly varying set points used in the cascade).15 This feature precludes the direct derivation of the underlying ODE representation; a way to overcome this problem involves the design of a “precompensator” that modifies the original DAE model into a new DAE system whose state space is control-independent. In ref 16, such a dynamic state feedback precompensator was proposed. However, this approach requires access to the state measurements, which clearly limits its applicability. In the present work, we establish that, for integrated process networks, inherent rank properties of the DAEs modeling the slow dynamics allow the derivation of a dynamic precompensator that does not require state measurements and leads to a DAE model with a controlindependent state space. Specifically, the precompensator consists of a minimal number of integrators added to suitable manipulated input channels. An ODE representation of the modified DAE model can then be derived, and an output feedback supervisory controller can be designed on the basis of this ODE model. The method is applied to two chemical engineering examples. The first is a gas-phase reactor connected to a condenser through large material recycle. The second application concerns a high-purity distillation column with a large recycle flow rate. Control of a high-purity distillation column is a challenging problem owing to highly nonlinear behavior, ill-conditioning, and a strong coupling between the top and bottom of the column.17 We focus, in particular, on a two-point control problem in such a column, i.e., control of both the bottom and top compositions (known to be especially challenging).
10.1021/ie0341909 CCC: $27.50 © 2004 American Chemical Society Published on Web 06/09/2004
Ind. Eng. Chem. Res., Vol. 43, No. 14, 2004 3529
1 x˘ ) f(x) + gs(x) us + gl(x) ul
Figure 1. Process network.
For both examples, simulations illustrate the performance and robustness of the proposed controllers. 2. Integrated Process Networks We consider the generic network shown in Figure 1, consisting of N processes (e.g., reactors, separation systems) and a material recycle stream for which the flow rate FR is much larger than the feed flow rate F0. Note that a large recycle flow rate FR implies equally large internal flow rates (F1, ..., FN-1). The presence of large internal and recycle flow rates compared to the external flow rates induces a time scale separation in the network dynamics, with the dynamics of the individual units evolving in a fast time scale and the dynamics of the overall network in a slow time scale.15 The mathematical model describing the overall and component material balances of such networks has the form N-1
x˘ ) f(x) +
gj(x) Fj + gR(x) FR ∑ gj(x) Fj + ∑ j)0,N j)1
where x is the vector of state variables (x ∈ X ⊂ Rn) and f(x), gj(x), and gR(x) are n-dimensional vector functions. Note that the two last terms involve large parameters corresponding to the internal flow rates F1, ..., FN-1 and the recycle flow rate FR. The effect of these large parameters can be isolated in a single parameter defined as the ratio FRs/F0s, where the subscript s denotes steady-state values. Extracting this ratio from these terms, we obtain
x˘ ) f(x) +
∑ gj(x) Fj +
j)0,N
FRs F0s
[
]
Fj FR gj(x) + gR(x) F0s FRs FRs j)1
N-1
F0s
∑
or, equivalently,
x˘ ) f(x) +
Fj
∑ gj(x) FjsF j)0,N FRs F0s
[∑
+
js
]
Fjs Fj FR gj(x) F0s + gR(x) F0s (1) FRs Fjs FRs j)1
N-1
When the ratios uj ) Fj/Fjs, for j ) 1, ..., N, R, and the small parameter ) F0s/FRs are defined, eq 1 takes the form
x˘ ) f(x) +
∑ gj(x) Fjsuj +
j)0,N
[∑
]
Fjs gj(x) F0s uj + gR(x) F0suR (2) j)1 FRs
1
N-1
where the stiffness is isolated in the small parameter and all of the other terms are of the same order of magnitude. Equation 2 can be written in a more compact form as
(3)
where us ∈ Rms is a vector of scaled input variables corresponding to the small flow rates (F0 and FN), ul ∈ Rml is a vector of scaled input variables corresponding to the large flow rates (FR and Fj for j ) 1, ..., N - 1), and gs(x) and gl(x) are matrices of appropriate dimensions. Note that the system in eq 3 is not in a standard singularly perturbed form because there is no explicit separation between fast and slow variables. In general, all of the variables x exhibit a fast initial transience followed by slow dynamics. In the fast time scale τ ) t/, in the limit f 0, the dynamics of the system take the form
dx/dτ ) gl(x) ul
(4)
Note that only the variables ul associated with the large flow rates are available for control purposes in this time scale. These typically include stabilization of holdups acting as integrators, which can be achieved by using simple proportional controllers
ul ) K(yl - ylsp) + I
(5)
where ylsp denotes the set point for the outputs yl ∈ Rml (the holdups that need to be stabilized), K is a ml × ml diagonal matrix, and I is a ml-dimensional vector whose entries are the unity I ) [1, ..., 1]T. Then, the fast dynamics of the network takes the form
dx/dτ ) gl(x) [K(yl - ylsp) + I]
(6)
The system in eq 3 is then written as
1 x˘ ) f(x) + gs(x) us + gl(x) [K(yl - ylsp) + I]
(7)
In the slow time scale t, when eq 7 is multiplied by and considering the limit f 0, the following quasisteady-state constraints are obtained
0 ) gl(x) [K(yl - ylsp) + I]
(8)
Note that these constraints are not linearly independent.15 In general, for a recycle network with C components in each process unit, n - C linearly independent terms can be isolated, leading to the following reformulation of the constraints in eq 8
0 ) gl(x) [K(yl - ylsp) + I] ) b(x) g j l(x) [K(yl - ylsp) + I]
(9)
where b(x) and g j l(x) are matrices of dimensions n × p and p × ml, where p ) n - C, b(x) is full column rank, and g j l(x) is full column rank. Note that the components of the vector g j l(x) [K(yl - ylsp) + I] are linearly independent by construction, which implies that the Jacobian of that vector is full row rank. Defining z ) limf0 {g j l(x) [K(yl - ylsp) + I]}/, the slow dynamics of the network take the form
x˘ ) f(x) + gs(x) us + b(x) z 0)g j l(x) [K(yl - ylsp) + I]
(10)
3530 Ind. Eng. Chem. Res., Vol. 43, No. 14, 2004
Equation 10 is a high-index DAE system because the solution for the (algebraic) variables z cannot be obtained directly from the algebraic equations. Let ys denote the output variables that are associated with control objectives in the slow time scale (e.g., product composition, total holdup). Note that, in the slow time scale, the small flow rates us are available as manipulated inputs. Typically, some of the set points of the fast output variables, ylsp, are also used as additional manipulated inputs. Such a cascaded control configuration becomes necessary when the number of controlled outputs ys in the slow time scale exceeds the number of available variables in us. In the case of such configurations, the algebraic constraints in the DAE description of the slow dynamics explicitly involve the manipulated input variables, which leads to a DAE model with control-dependent state space. Without loss of generality, it can be assumed that the m first components of the vector ylsp correspond to the set points used as additional manipulated inputs. After splitting the matrix K ) [K1 K2] where K1 is a ml × m matrix, the constraints of the DAE model of eq 10 yield l 0)g j l(x) [Kyl - K2ysp + I] - g j l(x) K1usp 2
(11)
l where usp ) [ysp ] is the m-dimensional vector of set 1 points used as additional manipulated inputs. Note that the matrix K1, as a submatrix of the diagonal matrix K, is a full column rank matrix where each column has only one constant nonzero entry. This feature, coupled with the matrix g j l(x) being full column rank, implies that the matrix g j l(x) K1 is also full column rank. In addition to being full column rank, the matrix g j l(x) K1 has m1 constant rows that correspond to the linearly independent constraints originating from the overall material balances. Note that an upper bound for m1 is the number of outputs that need to be stabilized in the fast time scale (e.g., holdups), i.e., ml. For the process networks under consideration, the following property is also assumed to hold: the p × m matrix g j l(x) K1 has a m × m constant and invertible submatrix (m e m1 e ml). Note that this assumption is trivially satisfied when only one set point is used as an additional manipulated input (i.e., m ) 1), which is often the case in practical applications (see the case studies in the present paper). This assumption implies that eq 11 can be multiplied by an invertible matrix such that the constraints of the DAE model of eq 10 take the form
where the following rank properties are fulfilled: (i) the n × p matrix b(x) is full column rank, (ii) C ˆ is a constant and invertible m × m (with m < p) matrix that accounts for the use of m of the set points ylsp as additional manipulated inputs, and (iii) the Jacobian of the vector [k h (x)T, k(x)T]T is full row rank. In what follows, we will consider DAE systems in the form of eq 12. The dependence of the state space on the manipulated inputs precludes a direct derivation of an expression for the algebraic variables by differentiation of the constraints. However, the rank conditions mentioned earlier allow us to modify the DAE system into a new DAE system of index 2 with a control-independent state space. 3. Dynamic Precompensator For the DAE system of the form in eq 12, we seek a dynamic precompensator with the goal of modifying the constraints that involve the m set points used as additional manipulated inputs in order to obtain a modified DAE system for which (i) the resulting constraints are independent of the new manipulated inputs and (ii) upon one differentiation the new constraints are solvable in the algebraic variables z; i.e., the modified DAE system has an index of 2. The above characteristics imply that no additional constraints will be present and, hence, the modified DAE system will have a controlindependent state space. The following proposition summarizes the main theoretical result of the paper. Proposition. Consider the DAE system
x˘ ) f(x) + b(x) z + gs(x) us 0)k h (x) + C ˆ usp 0 ) k(x)
(13)
where x ∈ Rn is the vector of state variables, z ∈ Rp the vector of algebraic variables, [us, usp]T ∈ Rms+m the vector of manipulated inputs, ys ∈ Rms+m the vector of controlled outputs, f(x), k h (x), and k(x) are n-dimensional, mdimensional, and (p - m)-dimensional vectors, respecˆ are matrices of appropriate tively, and b(x), gs(x), and C dimensions. Assuming also that (i) the n × p matrix b(x) is full column rank, (ii) C ˆ is a constant and invertible m × m (with m < p) matrix, and (iii) the Jacobian of the vector [k h (x)T, k(x)T]T is full row rank, then the dynamic extension
0)k h (x) + C ˆ usp
w˘ ) v1
0 ) k(x)
us ) v2
where C ˆ is a constant invertible m × m matrix (with m < p). Moreover, given that the Jacobian of the vector l + I] has full row rank, the Jacobian g j l(x) [Kyl - K2ysp 2 of the vector [k h (x)T, k(x)T]T is also full row rank. Therefore, the slow dynamics of the process networks under consideration are modeled by the following DAE system:
where usp ) w yields an index-2 DAE system whose state space is independent of the vector of manipulated inputs v ) [v1T, v2T]T. Proof. The direct substitution of the dynamic extension of eq 14 into the DAE system of eq 13 yields the following DAE system:
s
x˘ ) f(x) + b(x) z + g (x) u
]
f(x) b(x) 0 gs(x) x˘ ) + z+ v Im 0 w˘ 0 0
s
0)k h (x) + C ˆ usp 0 ) k(x)
[][ ][ ] [
(14)
0)k h (x) + C ˆw (12)
0 ) k(x)
(15)
Ind. Eng. Chem. Res., Vol. 43, No. 14, 2004 3531
The state space of this modified DAE system is clearly independent of the manipulated inputs v. Upon one differentiation of the constraints, the matrix coefficient for the algebraic variables takes the following form:
B)
[] ] ∂k h ∂x b(x) ∂k ∂x
(16)
Given the properties i and iii, it is clear that the matrix B as the product of a full row rank matrix and a full column rank matrix is invertible; the constraints are now solvable in z:
z ) -[Lb(x)k(x)]-1 [Lf(x)k(x) + Lgs(x)k(x) v2 + C h v1] ) R(x) + S1(x) v1 + S2(x) v2
(17)
where
k(x) )
[ ]
[]
k h (x) C ˆ , C h ) 0 k(x)
Lb(x)k(x) )
[] ] ∂k h ∂x ∂k b(x) ∂x
xj˘ ) hf (xj) + g j (xj) v
(18)
is a state space realization of the modified DAE system, where xj ) [xT, wT]T, the extended state vector, is constrained to evolve on the manifold defined by the constraints in eq 15, v is the new manipulated input vector, and
[
]
f(x) + b(x) R(x) , 0 b(x) S1(x) b(x) S2(x) + gs(x) g j (xj) ) Im 0
[
This is in contrast to the previous result in ref 16, where a dynamic state feedback compensator was proposed as a means of dealing with the dependence of constraints on the manipulated inputs. Remark 3. Note that, because the matrix g j l(x) K1 is full row rank, it admits a m × m invertible submatrix. For process networks for which the invertible submatrix is not constant, the constraints in eq 11 can be multiplied by an invertible matrix such that the slow dynamics of such process networks take the form
x˘ ) f(x) + b(x) z + gs(x) us
This implies that the DAE system in eq 15 is of index 2, which completes the proof. A direct substitution of the solution for z in the differential equation for x yields a state space realization of the original system. This is given in the following corollary. Corollary. Consider a DAE system of the form in eq 13 for which the rank conditions i-iii are satisfied, subject to the dynamic precompensator of eq 14. Then the dynamic system
hf (xj) )
Figure 2. Reactor-condenser network.
]
(19)
where z ) R(x) + S1(x) v1 + S2(x) v2. Remark 1. Note that the state variable xj is constrained to evolve on the p-dimensional manifold defined by the algebraic equations in eq 15. Therefore, the slow dynamics of the network is of order n + m - p, while the fast dynamics is of dimension p. Remark 2. Note that the proposed compensator consists of a dynamic extension whereby a minimal number of integrators are added, more specifically only to the manipulated input channels corresponding to the m set points used as additional manipulated inputs.
0)k h (x) + C ˆ (x) usp 0 ) k(x) where the n × p matrix b(x) is full column rank, the m × m matrix C ˆ (x) is invertible, and the Jacobian of the vector [k h (x)T, k(x)T]T is full row rank. For such DAE systems, the derivation of a dynamic precompensator to obtain an index-2 modified DAE system depends on the structure of the matrix C ˆ (x) and on the outputs to be controlled.18 On the basis of the state space realization in the previous corollary, an output feedback controller can be designed using existing techniques for nonlinear ODE systems.19,20 4. Case Studies 4.1. Reactor-Condenser Network. We consider a gas-phase reactor, which is fed with pure reactant A at a constant molar flow rate F0. The following two firstorder reactions yield the high boiling point desired product C: k1
k2
A 98 B 98 C The reactor effluent F is fed to a separation unit, where the product C condenses. The product-rich liquid phase of the condenser is withdrawn from the condenser bottom at a flow rate L, while the gas phase from the condenser is recycled to the reactor at a flow rate R (see Figure 2). We assume that, to enhance conversion, the recycle flow rate is significantly larger than the feed flow rate.
3532 Ind. Eng. Chem. Res., Vol. 43, No. 14, 2004
A dynamic model of this network is given by the following ODE system: M ˙ R ) F0 + R - F 1 y˘ A,R ) [F (y - yA,R) + R(yA,C - yA,R) - k1MRyA,R] MR 0 A,0 1 y˘ C,R ) [-F0yC,R + R(yC,C - yC,R) + k2MR(1 - yA,R - yC,R)] MR
M ˙ V)F-R-N 1 [F(yA,R - yA,C) - NA + yA,CN] y˘ A,C ) MV 1 [F(yC,R - yC,C) - NC + yC,CN] y˘ C,C ) MV
M ˙ L)N-L 1 [N - xAN] x˘ A ) ML A 1 [N - xCN] x˘ C ) ML C
}
}
}
reactor
condenser (vapor phase)
condenser (liquid phase)
where yA,R and yC,R are the mole fractions of components A and C in the reactor, yA,C and yC,C the mole fractions of components A and C in the condenser, and xA and xC are the corresponding liquid mole fractions in the condenser. The interphase mass-transfer flow rate N is expressed as the sum of the mass flow rates of the components (i.e., N ) NA + NB + NC) with
[
Nj ) KjA yj,C -
pjs
]
xj
ptot
ML FL
[ ]
R h - κF h R h (y - yA,R) MR A,C R h (yC,C - yC,R) M gl(x) ul ) F0 R κF h -R h κF h (y - yA,C) MV A,R κF h (y - yC,C) MV C,R
Note that, in this time scale, the compositions and the holdup in the condenser liquid phase remain constant because the large flow rates only affect the reactor and condenser gas-phase holdups and compositions. Thus, in this time scale, only the gas holdups MR and MV need to be stabilized, which is easily achieved by using simple proportional controllers involving the outputs ul:
R h )1-K h p1(MVs - MV) F h )1-K h p2(MRs - MR)
j ) A-C
where Kj is the mass-transfer coefficient for component j expressed per unit of interfacial area, FL the molar density of the condenser liquid phase, A the interfacial area per unit liquid holdup volume, ptot the total pressure in the condenser unit, and pjs the saturation vapor pressure for component j at the temperature of the condenser Tcond. The Antoine equation is used to evaluate the saturation vapor pressure
pjs ) eaj-bj/(Tcond+cj)
subscript s refers to the nominal steady-state values and O(‚) is the standard order of magnitude notation, are defined, the process model takes the general form of eq 7, where x is the vector of state variables (compositions and holdups in the reactor and condenser), us ) [L] ∈ R is the manipulated input corresponding to the small flow rate (because the feed flow rate is assumed to be h, F h ]T ∈ R2 is the vector of constant), and ul ) [R manipulated inputs corresponding to the large flow h ) F/Fs. rates with R h ) R/Rs and F In the fast time scale (τ ) t/), in the limit f 0, the dynamics take the form of eq 4 where
(21)
In the slow time scale, the fast dynamics give rise to the constraints 0 ) gl(x) ul. The last three constraints corresponding to the condenser vapor phase are linear combinations of the other constraints. Defining
j ) A-C
where aj, bj, and cj are the Antoine constants for component j. The total pressure is obtained from the following algebraic equation:
Vcond ) VL + VV ) FLML +
MVRgasTcond ptot
(20)
where Vcond is the condenser volume, VL and VV are the volumes of the condenser liquid and vapor phases (respectively), and Rgas is the ideal gas law constant. The key outputs to be controlled are the holdups MR, MV, and ML (for the reactor and condenser, respectively), which behave like integrators, and the product purity xC. When the singular perturbation parameter ) F0/Rs and the parameter κ ) Fs/Rs ) O(1), where the
we have gl(x) ul ) b(x) g j l(x) ul, and the dynamics take the general form of eq 10. In this time scale, the outputs to be controlled are the composition of the desired product xC and the total holdup (MR + MV + ML), which
Ind. Eng. Chem. Res., Vol. 43, No. 14, 2004 3533
is not influenced by the large flow rates and hence does not evolve in the fast time scale. The small product flow rate L is the only available manipulated input, hence the necessity to introduce another manipulated input. To this end, a cascaded configuration was employed, where the set point for the reactor holdup MRs was used as an additional input. Considering L and MRs as the manipulated inputs, the DAE system for the slow dynamics can be expressed in the form of eq 12:
x˘ ) f(x) + b(x) z + gs(x) us 0)k h (x) + C ˆ usp 0 ) k(x) where us ) L is the manipulated input corresponding to the small flow rate and usp ) MRs is the set point used as a manipulated input
C ˆ ) [κ K h p 2] k h (x) ) {[1 - K h p1(MVs - MV)] - κ(1 + K h p2MR)} and the two linearly independent constraints that do not involve the inputs are
k(x) )
[
[1 - K h p1(MVs - MV)](yA,C - yA,R)
[1 - K h p1(MVs - MV)](yC,C - yC,R)
]
Note that the constraint involving the manipulated input, i.e., k h (x) + C ˆ usp, relates the reactor holdup MR and the condenser vapor holdup MV. Moreover, eq 20 introduces a relationship between the condenser holdups, MV and ML. Thus, controlling the condenser liquid holdup ML suffices to stabilize the total holdup. Therefore, the controlled outputs are y1 ) ML and y2 ) xC. It is clear that the rank conditions i-iii are fulfilled, which leads to the following precompensator of the form in eq 14:
w˘ ) v1 us ) L ) v2 usp ) MRs ) w The dynamic extension corresponds to adding an integrator to the channel of the manipulated input MRs. The resulting DAE system leads to a state space realization that is unstable because of the condenser liquid holdup behaving as an integrator. The condenser liquid holdup is first stabilized by means of a proportional controller involving the second manipulated input, i.e., the small flow rate L. An output feedback controller with integral action is then designed to control the product composition. The relative order of y2 with respect to v1 is r ) 2 so that a controller is designed to enforce a critically damped second-order response
dy2 d2y2 y2 + γ21 + γ22 2 ) y2,sp dt dt where y2,sp denotes the set points for the product composition. The controller consists of an input/output
Table 1 variable F0 R F L MR MV ML k1 k2 Kp1 Kp2 Kp3 Tcond Treactor Vcond Vreactor A Rgas FL KA KB KC aA, bA, cA aB, bB, cB aC, bC, cC
description feed flow rate (mol/s) recycle flow rate (mol/s) recycle flow rate (mol/s) withdrawing flow rate (mol/s) reactor holdup (mol) condenser vapor holdup (mol) condenser liquid holdup (mol) reaction rate (s-1) reaction rate (s-1) proportional controller gain (s-1) proportional controller gain (s-1) proportional controller gain (s-1) condenser temperature (K) reactor temperature (K) condenser volume (m3) reactor volume (m3) interphase area (m2/m3) ideal gas law constant (J/mol‚K) condenser liquid density (mol/m3) mass-transfer coefficient for A (mol/m2‚s) mass-transfer coefficient for B (mol/m2‚s) mass-transfer coefficient for C (mol/m2‚s) Antoine coefficients a, b, and c for component A Antoine coefficients a, b, and c for component B Antoine coefficients a, b, and c for component C
value 112 1290 1402 112 5215 436 6924 0.06 0.029 15 15 0.5 279 373 6.00 5.69 100 8.314 34 54 889 30 30 16 27.3, 4400, 4 27.4, 4100, 0 19.1, 3400, 0
linearizing state feedback controller coupled with an “open-loop” observer and an external linear controller.21 The controller was tuned with the parameters γ21 ) 600 s and γ22 ) 90 000 s2. Table 1 gives the nominal values of the process variables at steady state. Figure 3 shows the performance of this controller, for a 5.45% increase in the purity set point of component C, xC, at t ) 20 min. The controller exhibits a good tracking performance, imposing the requested critically damped second-order response. Figure 4 shows the performance of the controller for the same set-point change, in the presence of -5% error in the masstransfer coefficient KC and -5% error in the reaction rate k2. The controller yields a good performance in tracking the set-point change and rejecting the effects of these errors. 4.2. High-Purity Distillation Column. We consider a distillation column with N trays (numbered from top to bottom), to which a saturated liquid containing a mixture of three components with mole fractions x1f and x2f of components 1 and 2, respectively, is fed at (molar) flow rate F on tray Nf. The heavy component 3, which is the desired product, is removed at the bottom from the reboiler at a flow rate B, while the lighter components 1 and 2 are removed at the top from the total condenser at a flow rate D. In this column, a large vapor boilup VB and liquid recycle R are used compared to the feed, distillate, and bottom product flow rates to attain a high purity of the desired component 3 in the bottom product. Constant molar overflow and constant holdup on each tray (at possibly different values) are assumed. Under the above assumptions, a standard dynamic model of the column is obtained, which is given by the following ODE system:22
3534 Ind. Eng. Chem. Res., Vol. 43, No. 14, 2004
Figure 3. Reactor-condenser network: closed-loop input/output responses, in the nominal case, for a 5.45% increase in the purity set point of component C, at t ) 20 min.
M ˙ C ) VB - R - D VB (y - x1,D) x˘ 1,D ) MC 1,1 VB x˘ 3,D ) (y - x3,D) MC 3,1
}
condenser
1 [V (y - y1,i) + R(x1,i-1 - x1,i)] Mi B 1,i+1 1 x˘ 3,i ) [VB(y3,i+1 - y3,i) + R(x3,i-1 - x3,i)] Mi
x˘ 1,i )
}
tray i < Nf
1 [V (y - y1,i) + R(x1,i-1 - x1,i) + F(x1f - x1,i)] Mi B 1,i+1 1 x˘ 3,i ) [VB(y3,i+1 - y3,i) + R(x3,i-1 - x3,i) + F(x3f - x3,i)] Mi
x˘ 1,i )
} }
feed tray i ) Nf 1 [V (y - y1,i) + R(x1,i-1 - x1,i) + F(x1,i-1 - x1,i)] Mi B 1,i+1 1 x˘ 3,i ) [VB(y3,i+1 - y3,i) + R(x3,i-1 - x3,i) + F(x3,i-1 - x3,i)] Mi
x˘ 1,i )
}
tray i > Nf M ˙ R ) R - VB + F - B 1 x˘ A,B ) [R(x1,N - x1,B) - VB(y1,B - x1,B) + F(x1,N - x1,B)] MR 1 x˘ 3,B ) [R(x3,N - x3,B) - VB(y3,B - x3,B) + F(x3,N - x3,B)] MR
reboiler where MC, x1,D, x3,D, y1,D, and y3,D are the molar liquid holdup, liquid mole fractions, and vapor mole fractions of components 1 and 3 in the condenser, Mi, x1,i, x3,i, y1,i, and y3,i are the molar liquid holdup, liquid mole
fractions, and vapor mole fractions of components 1 and 3 in tray i, and MR, x1,B, x3,B, y1,B, and y3,B are the corresponding holdup, liquid mole fractions, and vapor mole fractions in the reboiler. The key outputs to be controlled are the bottom purity x3,B, the top composition x1,D, and the two liquid holdups MC and MR that behave like integrators. When the singular perturbation parameter ) Ds/Rs and κ ) VBs/ Rs ) O(1) are defined, the process model takes the general form of eq 7,22 where x is the vector of state variables (compositions and holdups in each stage), us ) [D, B]T ∈ R2 is the vector of manipulated inputs corresponding to small flow rates, and ul ) [R h, V h B] T ∈ R2 is the vector of manipulated inputs corresponding to large flow rates where R h ) R/Rs and V h B ) VB/VBs. In the fast time scale (τ ) t/), in the limit f 0, the dynamics take the form of eq 4, where
Ind. Eng. Chem. Res., Vol. 43, No. 14, 2004 3535
[
]
Figure 4. Reactor-condenser network: closed-loop input/output responses for a 5.45% increase in the purity set point of component C, in the presence of -5% error in the mass-transfer coefficient KC and -5% error in the reaction rate k2, at t ) 20 min.
for 1 e i e N. In this time scale, only the outputs ul are available for control purposes. In particular, the control of the liquid holdups in the condenser and reboiler (MC and MR) is easily achieved by using simple proportional controllers:
R h )1-K h c1(MCs - MC) V hB ) 1 - K h c2(MRs - MR) In the slow time scale t, the fast dynamics give rise to the constraints 0 ) gl(x) ul. The last three constraints (corresponding to the reboiler) are linear combinations of the other ones. Defining
-1 0 0 l g j (x) ) x1,i-1 - x1,i x3,i-1 - x3,i l
κ κ(y1,1 - x1,D) κ(y3,1 - x3,D) κ(y1,i+1 - y1,i) κ(y3,i+1 - y3,i) l
for 1 e i e N, we have gl(x) ul ) b(x) g j l(x) ul, and the dynamics take the general form of eq 10. Only the small flow rates D and B affect the slow dynamics. In this time scale, the outputs to be controlled consist of the top and bottom compositions, x1,D and x3,B, respectively, and the total liquid holdup (MC + MR), which needs to be stabilized because it is not affected by the large flow rates.22 Hence, we need an additional manipulated input. A natural approach to this end is a cascaded control configuration where one of the set points for the condenser-reboiler holdups used in the fast proportional control is treated as an additional manipulated input variable. Considering MCs as an additional manipulated input, the DAE system for the slow dynamics can be expressed in the form of eq 12:
x˘ ) f(x) + b(x) z + gs(x) us 0)k h (x) + C ˆ usp 0 ) k(x) where us ) [D, B]T is the vector of manipulated inputs corresponding to the small flow rates and usp ) [MCs] is the set point used as manipulated input
C ˆ ) [K h c1] h c2(MRs - MR)]} k h (x) ) {-(1 + K h c1MC) + κ [1 - K and
and the 2N + 2 linearly independent constraints that do not involve the inputs are
3536 Ind. Eng. Chem. Res., Vol. 43, No. 14, 2004
Figure 5. High-purity distillation column: closed-loop input/output responses, in the nominal case, for a 2.35% increase in the bottom purity x3,B and a 22.4% decrease in the top composition x1,D, at t ) 15 min.
Figure 6. High-purity distillation column: closed-loop input/output responses for a 2.35% increase in the bottom purity x3,B and a 22.4% decrease in the top composition x1,D, in the presence of 2% error in R1 and -2% and +3% step disturbances in the feed composition x1f and x2f, at t ) 15 min.
[
k(x) ) κ[1 - K h c2(MRs - MR)](y3,1 - x3,D)
κ[1 - K h c2(MRs - MR)](y1,1 - x1,D)
κ[1 - K h c2(MRs - MR)](y1,i+1 - y1,i + x1,i-1 - x1,i)
κ[1 - K h c2(MRs - MR)](y3,i+1 - y3,i + x3,i-1 - x3,i) l
]
where 1 e i e N, x1,i and x3,i are the liquid mole fractions of 1 and 3 in tray i, and y1,i and y3,i are the vapor mole fractions in tray i. It can easily be verified that the rank conditions necessary for the application of the proposition are fulfilled, so that the following precompensator of the form in eq 17 is obtained:
Ind. Eng. Chem. Res., Vol. 43, No. 14, 2004 3537 Table 2 variable
description
value
B D F Kc1 Kc2 Kc3 MC Mi MR N Nf R VB R1 R2
bottom product flow rate (mol/min) distillate flow rate (mol/min) feed flow rate (mol/min) proportional controller gain (min-1) proportional controller gain (min-1) proportional controller gain (min-1) condenser liquid holdup (mol) liquid holdup on tray i (mol) reboiler liquid holdup (mol) total number of trays feed tray liquid recycle flow rate (mol/min) vapor boilup flow rate (mol/min) relative volatility of component 1 relative volatility of component 2
50.0 50.0 100.0 20.0 20.0 5.0 180.0 175.0 200.0 15 8 1000.0 1050.0 1.5 1.3
w˘ ) v1 us )
[][] v D ) v2 B 3
usp ) MCs ) w
Acknowledgment
The dynamic extension corresponds to adding an integrator to the channel of the manipulated input MCs. The resulting DAE system leads to a state space realization, which is unstable owing to the total holdup that behaves as an integrator. The total holdup is first stabilized by means of a proportional controller involving the second manipulated input, i.e., the small flow rate D. An output feedback controller with integral action is designed to control the top and bottom compositions. The relative orders of these two outputs are r1 ) r2 ) 1 with respect to the new manipulated inputs v, so that a controller is designed to enforce the response
yi + γi
dyi ) yi,sp, i ) 1, 2 dt
control configuration, the slow dynamics are modeled by a high-index DAE whose state space is controldependent. By using inherent rank properties of such networks, we proposed a dynamic precompensator to obtain a new DAE model with a control-independent state space. Specifically, a minimal number of integrators is added to the manipulated input channels corresponding to the set points used in the cascade. The design of a controller addressing the objectives in the slow time scale was then performed in two steps. The first step includes the dynamic precompensation and subsequent derivation of an ODE representation of the resulting DAE model, and the second step involves the design of a nonlinear output feedback controller on the basis of the underlying ODE system. The procedure was applied to two chemical engineering examples, a reactor-condenser network, and a distillation column. In both examples, a shortage of manipulated inputs in the slow time scale necessitated a cascade control structure; the proposed method allowed one to successfully address the nonlinear output feedback control of those systems.
(22)
where y1,sp and y2,sp denote the set points for the top and bottom compositions, respectively. The controller consists of an input/output linearizing state feedback controller coupled with an open-loop observer and an external linear controller.21 The controller was tuned with the parameters γ1 ) γ2 ) 20 min. Table 2 gives the nominal values of the process variables at steady state. Figure 5 shows the performance of this nonlinear output feedback controller, for a 2.35% increase in the bottom purity x3,B and a 22.4% decrease in the top composition x1,D at t ) 15 min. Note that the increase in the bottom purity and the decrease in the top composition have to lead to a consistent steady state. The controller performs excellent tracking. Figure 6 shows the performance of the controller for the same set-point changes, in the presence of 2% error in R1, -2% and 3% unmeasured step disturbances in the feed composition x1f and x2f, at t ) 15 min. Clearly, the controller eliminates the effect of disturbances and induces the desired input/output behavior. 5. Conclusion In this work, we considered integrated process networks with large recycle. The coexistence of small and large flow rates is known to induce a time scale multiplicity in the network dynamics. In a cascaded
Financial support by ACS-PRF and by NSF/CTS Grant 0234440 is gratefully acknowledged. Literature Cited (1) Grossman, I. E.; Westerberg, A. W. Research challenges in process systems engineering. AIChE J. 2000, 46, 1700-1703. (2) Marquardt, W. Fundamental modeling and model reduction for optimization based control of transient processes. Preprints of Chemical Process Control-6, Tucson, AZ, 2000; p 30. (3) Bildea, C. S.; Dimian, A. C.; Iedema, P. D. Nonlinear behavior of reactor-separator-recycle systems. Comput. Chem. Eng. 2000, 24, 209-214. (4) Pushpavanam, S.; Skogestad, S. Nonlinear behavior of an ideal reactor separator network with mass recycle. Chem. Eng. Sci. 2001, 57, 2837-2849. (5) Jacobsen, E. On the dynamics of integrated plants: nonminimum phase behavior. J. Process Control 1999, 9, 439-451. (6) Morud, J.; Skogestad, S. Analysis of instability in an industrial ammonia reactor. AIChE J. 1998, 44, 888-895. (7) Luyben, M. L.; Floudas, C. A. Analyzing the interaction of design and controls2. reactor-separator-recycle system. Comput. Chem. Eng. 1994, 18, 971-993. (8) Yi, C. K.; Luyben, W. L. Design and control of coupled reactor/column systemssparts 1-3. Comput. Chem. Eng. 1997, 21, 25-68. (9) Ng, C.; Stephanopoulos, G. Synthesis of control systems for chemical plants. Comput. Chem. Eng. 1996, 20, S999-S1004. (10) Luyben, M. L.; Tyreus, B. D.; Luyben, W. L. Plantwide control design procedure. AIChE J. 1997, 43, 3161-3174. (11) Wang, P.; McAvoy, T. Synthesis of plantwide control systems using dynamic model and optimization. Ind. Eng. Chem. Res. 2001, 40, 5732-5742. (12) Larsson, L.; Govatsmark, M. S.; Skogestad, S.; Yu, C. C. Control structure selection for reactor, separator and recycle processes. Ind. Eng. Chem. Res. 2003, 42, 1225-1234. (13) Farschman, C. A.; Viswanath, K. P.; Ydstie, B. E. Process systems and inventory control. AIChE J. 1998, 44, 1841-1857. (14) Hangos, K. M.; Alonso, A. A.; Perkins, J. D.; Ydstie, B. E. Thermodynamical approach to the structural stability of process plants. AIChE J. 1999, 45, 802-816. (15) Kumar, A.; Daoutidis, P. Nonlinear dynamics and control of process systems with recycle. J. Process Control 2002, 12, 475484. (16) Kumar, A.; Daoutidis, P. Control of non-linear differential algebraic equation systems; Chapman and Hall/CRC: London, 1999. (17) Skogestad, S. Dynamics and control of distillation columns, a critical survey. Model., Identif. Control 1997, 18, 177-217.
3538 Ind. Eng. Chem. Res., Vol. 43, No. 14, 2004 (18) Contou-Carrere, M.-N.; Daoutidis, P. Output feedback regularization of non-linear DAE systems. Proceedings of the 15th IFAC World Congress, Barcelona, Spain, 2002. (19) Kravaris, C.; Kantor, J. Geometric methods for nonlinear process control. parts 1-2. Ind. Eng. Chem. Res. 1990, 29, 22952323. (20) Isidori, A. Nonlinear control systems; Springer-Verlag: London, 1995. (21) Daoutidis, P.; Kravaris, C. Dynamic output feedback control of minimum-phase multivariable nonlinear processes. Chem. Eng. Sci. 1994, 49, 433-447.
(22) Kumar, A.; Daoutidis, P. Nonlinear Model Reduction and Control for High-Purity Distillation Columns. Ind. Eng. Chem. Res. 2003, 42, 4495-4505.
Received for review October 16, 2003 Revised manuscript received May 3, 2004 Accepted May 4, 2004 IE0341909