- Email: [email protected]
Self-oscillating and chaotic behaviour of a PI-controlled CSTR with control valve saturation
Journal of Process Control 14 (2004) 51–59 www.elsevier.com/locate/jprocont
Self-oscillating and chaotic behaviour of a PI-controlled CSTR with control valve saturation Manuel Pe´reza,1, Pedro Albertosb,* a
Department of ‘‘Fı´sica, Ingenierı´a de Sistemas y Teorı´a de la Sen˜al’’, Escuela Polite´cnica Superior, Alicante University, Campus de San Vicente, Alicante, Spain b Department of Systems Engineering and Control, E.T.S.I.I., Universidad Polite´cnica de Valencia, Camino de Vera s/n, Valencia, Spain Received 21 May 2002; received in revised form 16 April 2003; accepted 26 April 2003
Abstract The time response of a PI-controlled continuous stirred tank reactor (CSTR), where a first order irreversible reaction occurs, is analysed in this paper. The parameters of the PI controller and the saturation of the cooling water temperature flow are used as parameters to show that self-oscillations and chaotic dynamics can appear. The conditions under which self-oscillation occurs are investigated taking into account heat generation, heat removal and the variation of eigenvalues of the linearized system at the equilibrium points, respectively. Also, the appearance of a Silnikov type homoclinic orbit, giving chaotic dynamics, is realised. The existence of a new set of strange attractors with PI control has been verified and the sensitivity to the initial conditions is also analysed. # 2003 Elsevier Ltd. All rights reserved. Keywords: Control of continuous reactor; Steady state; Self-oscillation; Nonlinear and chaotic behaviour
1. Introduction The study of a continuous stirred tank reactor CSTR where a simple irreversible reaction A ! B occurs is basic for the analysis and control of chemical reactors. In the early years, this simple reaction scheme was considered as an ‘‘academic work’’ by many people, especially from the industry. However, the dynamics and the control of many reactors can be treated as processes in which the CSTR model is applicable. Processes involving catalytic reactors and hydrolysis reactions are, among others, good examples. The issues about the dynamic behaviour of CSTRs have been treated in various papers, [1–3], and books, and the oscillatory dynamics of a stirred tank reactor with a single reaction has been investigated using new numerical procedures in Kevrekidis [4]. In the past, various researchers [5,6] investigated a single irreversible * Corresponding author. Tel.: +34-96-387-9570; fax: +34-96-3879579. E-mail addresses: [email protected] (M. Pe´rez), [email protected] (P. Albertos). 1 Fax: +34-96-590-3682. 0959-1524/03/$ - see front matter # 2003 Elsevier Ltd. All rights reserved. doi:10.1016/S0959-1524(03)00032-5
exothermic reaction in which the coolant temperature was varied, following sine and square waves. They showed that with the variation of a single parameter, i.e. the amplitude of the variation in the coolant temperature, the reactor changes from quasi-periodic to periodic behaviour, and through a sequence of period doubling, it bifurcates to chaos. No control systems were used in any of these previous research studies. This is probably due to the fact that non-isothermal CSTR can show three steady states, not always being stable. The reactor evolves to one of its stable steady states [2] and without control, it is impossible to reach any arbitrary desired set point. Forced non-linear oscillators have been considered in various contexts. For example, the autonomous van der Pol oscillator, which is taken from electronics, is treated in Cartwright [7] showing the appearance of chaotic behaviour. Chaos can also occur via periodic perturbations in the Duffing equation [8]. Other examples can be found in Lichtenberg [9] and Seydel [10]. The simplest proportional controllers have been used to control the tank level and the temperature by manipulating the flow of the liquid leaving the tank and the cooling water [11–14]. In the paper of Pellegrini [15],
52
M. Pe´rez, P. Albertos / Journal of Process Control 14 (2004) 51–59
in order to get the chaotic behaviour of a PI-controlled CSTR, the outlet temperature is assumed as the manipulated variable without considering the cooling process dynamics. Moreover, more attention is paid to the non linear dynamics of the CSTR whereas the effect of the controller is not fully investigated. Also, in these conditions, it seems to be easier to reach the chaotic behaviour [6]. More recently, other more sophisticated control strategies for polymerization reactors based on feedback linearization and non-linear predictive control [16], or state observers [3] have been reported. In the work of Alvarez-Ramirez [17,18] the problem of the robust stabilisation of a CSTR has been studied and Femat [19], using robust asymptotic feedback, has corroborated the existence of the chaotic behaviour. In other papers, the self-oscillation in a CSTR using the flow rate as manipulated variable has been studied [20]. The PI control strategy for a CSTR provides a good example to analyse the difficulties encountered in controlling a process with several steady states. In this paper, an exothermic first order irreversible reaction in a CSTR is considered. The control system is very simple and common in the industry being composed by two PI controllers. The first controller manipulates the flow rate of the liquid leaving the reactor as a function of the volume in the tank error. It does not influence the qualitative response of the reactor and it can be assumed that it works perfectly, keeping constant the reactor volume. A second PI controller manipulates the flow rate of cooling water to the jacket as a function of the error in the reactor’s temperature. The paper shows that self-oscillation and chaotic behaviours are possible without any external periodic excitation. Just with proportional feedback control the reactor may reach self-oscillations independently of the values of the proportional gain. Another typical situation, like the saturation in the refrigerating control valve, is considered. This saturation may lead to a chaotic behaviour. Simulation studies are used for corroborating the existence of a homoclinic orbit of a special type, the so called Silnikov orbit [21], which contains an infinite number of unstable periodic orbits, that is, a chaotic behaviour.
2. CSTR model equations Let us assume an exothermic reaction A ! B in a CSTR. The reaction heat is removed by a cooling jacket surrounding the reactor. A first-order reaction in reactant A, negligible heat losses and perfectly mixing are assumed [11–13,22,23]. The jacket water is assumed to be perfectly mixed and the mass of the metal walls is considered negligible, so the thermal inertia of the metal is not considered. The reactor model without control includes the following equations:
dV ¼ Fo F dt
ð1Þ
d ðV C a Þ ¼ Fo Cao F Ca V Ca eE=RT dt
ð2Þ
dðVhÞ ¼ ðFo ho FhÞ HVCa eE=RT UA T Tj dt h ¼ cp T ; ho ¼ cp To ð3Þ
j Vj
dhj ¼ j Fj hjo hj þ UA T Tj dt hj ¼ cpj Tj ; hjo ¼ cpj Tjo
ð4Þ
obtained by a total mass balance (1), a component A mass balance (2), an energy balance in the reactor (3), and an energy balance in the jacket (4). A description of the variables and parameter values is given in Table 1. In addition to Eqs. (1)–(4), there are two PI controllers for the reactor outflow, F, and the cooling water flow, Fj, as shown in Fig. 1. The first one is acting on the inlet reactor flow to keep the volume constant. Note that the reactor volume, V(t), via manipulation of the outlet stream, is not affected by the rest of variables. Thus, in order to simplify the analysis, a perfect control of this volume is assumed, keeping it constant and equal to Vs. On the other hand, the volume does affect the temperature, the concentration and, thus, the refrigerating equations, and could be considered as a disturbance. The second one has the following equation: dFj dT 1 þ ðTr TÞ ¼ Kp ð5Þ dt Ki dt where Tr is the set point (operating point) of the reactor temperature and Kp, Ki are the controller proportional and integral action. The manipulated variable Fj is constrained: 0 4 Fj 4 Fjmax
ð6Þ
where Fjmax is the maximum value of the cooling water flow through the control valve CV2 (Fig. 1). The set of Eqs. (1)–(5) determine a good reactor model [13]. The variables Fo, To, Cao, and Tjo can be considered as load disturbances whose effect must be minimized by the control system. Remark 1. It will be shown that under nominal operating conditions, that is, V(t)=Vs and without any flow disturbance [Fo(t)=Fos], yielding VðtÞ=Fo ðtÞ constant, the reactor chaotic behaviour cannot be achieved despite other possible disturbances, such as To, Cao, and Tjo.
M. Pe´rez, P. Albertos / Journal of Process Control 14 (2004) 51–59
Remark 2. The dynamic related to the jacket temperature is much faster than that related to the reactor temperature. Thus, the jacket temperature time constant is negligible and a simplified model can be derived. With this assumption, from Eq. (4), it yields
Tj ¼
53
j cpj Fj Tjo þ UA T j cpj Fj þ UA
ð7Þ
To get the reduced model, Tj should be replace by this expression in Eq. (3).
Table 1 Variables, nominal operating conditions and parameter values Variable
Description
Value
V F Fo Ca h ho hj hjo T Tj Fj Fos Vs Cao Caor To Vj a E U A Tjo R H cp cpj r rj Tr
Reactor volume (m3) Volumetric flow rate for the outlet stream (m3/h) Volumetric flow rate for the inlet stream (m3/h) Reactant concentration for the outlet stream (kmol A/m3) Enthalpy for the outlet liquid (kJ/kg) Enthalpy for the inlet liquid (kJ/kg) Enthalpy for the outlet cooling water Enthalpy for the inlet cooling water Reactor temperature ( C) Jacket temperature ( C) Volumetric flow rate of cooling water (m3/h) Volumetric flow rate for the inlet stream (m3/h) Steady state reactor volume (m3) Reactant concentration inlet stream (kmol A/m3) Initial reactant concentration (kmol A/m3) Inlet stream temperature ( C) Jacket volume (m3) Preexponential factor from Arrhenius law (h1) Activation energy (kJ/kmol) Overall heat transfer in the jacket kJ/(h m2 C) Heat transfer area (m2) Inlet stream cooling water temperature ( C) Perfect-gas constant kJ/kmol K Enthalpy of reaction (kJ/kmol) Heat capacity inlet and out streams (kJ/kg K) Heat capacity of cooling water(kJ/kg K) Density of the inlet and out streams (kg/m3) Density of cooling water (kg/m3) Set point temperature ( C)
1.13 1.36 8 3.92 21.7 0.085 7.08 1010 69,815 3065 23.22 21.7 8.314 69,815 3.13 4.18 800 1000 60.6
Fig. 1. Perfectly mixed CSTR with two PI controllers.
54
M. Pe´rez, P. Albertos / Journal of Process Control 14 (2004) 51–59
2.1. CSTR dimensionless model To generalise the mathematical model of the reactor, the volume V being taken as a constant disturbance, Eqs. (2)–(5) can be expressed in a dimensionless way, leading to the following state-space model: dx2 x60 ¼ ðx20 x2 Þ c0 x2 e1=x3 d x 1 c2 dx3 x60 ¼ ðx30 x3 Þ þ c1 x2 e1=x3 ðx3 x4 Þ d x x 1 1 dx4 ¼ c x ð x x Þ þ c ð x x Þ 3 5 40 4 4 3 4 4 : d dx K x pd 5 60 ¼ ðx30 x3 Þ þ c1 x2 e1=x3 d x30 x1 c2 x3 x3r ðx3 x4 Þ þ x1 Kid
If the cooling water flow, Fjs, is not limited, providing as much refrigeration as required, the equilibrium points are obtained from (8), for x1e=1 and x3e ¼ x3r ;
x6e ¼ x60
ð10Þ
But in practice, as already noted in Eq. (6), this flow is limited and an equilibrium can be reached under saturated conditions. If this is the case, Fjs=Fjmax and the cooling water flow cannot decrease the reactor temperature, thus the equilibrium temperature reaches a value, Te, greater than the set point, (Te > Tr), and the corresponding reactant concentration will be smaller. The new equilibrium conditions are x5e ¼ x5max ;
x6e ¼ x6o
ð11Þ
Taking into account Eqs. (9) at the equilibrium points (11), the following relationship between the steady-state values can be obtained: ð8Þ
A full description of the dimensionless variables is given in Table 2. The fourth order model 4, (8), represents the reactor dynamics. If remark 2 is applicable, using Eq. (7) in (8), the reduced order model 3 will be dx2 x60 1=x3 d ¼ x1 ðx20 x2 Þ c0 x2 e dx3 x60 1=x3 d ¼ x1 ðx30 x3 Þ þ c1 x2 e c2 c3 x5 ðx3 x40 Þ 3 : x1 ðc3 x5 þ c4 Þ dx 5 ¼ KPd x60 ðx x Þ þ c x e1=x3 30 3 1 2 d x30 x1 c2 c3 x5 ðx3 x40 Þ þ ðx3 x3r Þ x1 ðc3 x5 þ c4 Þ Kid ð9Þ
3. Analysis of the effect of the control valve saturation It is well known that an exothermic CSTR without control has multiple steady states [24–27]. This is due to the fact that the curve of generated heat and the line of removed heat versus temperature can intersect at three points. A detailed behaviour of a CSTR without control can be found in Luyben [11], Stephanopoulos [12], Marlin [13] and Aris [25]. To run the reactor at any desired temperature T, the heat removal curve or the heat generation lines must be modified by changing some system parameters or adding a feedback control system. If a PI controller is used, Eqs. (2)–(4) should be considered but adding Eq. (5).
c2 c3 x5max x6e ðx30 x3e Þ x3e x40 c3 x5max þ c4 x6e x20 c1 ¼0 þ ðx3e x40 Þ ðx6e e1=x3e þ c0 Þ
ð12Þ
Eq. (12) shows the variation of the equilibrium dimensionless temperature as a function of the maximum value of dimensionless cooling water flow x5max. The plot of the implicit function f(x3e,x5max)=0, Eq. (12), is shown in Fig. 2. This plot shows the possibility of three steady states corresponding to the points P1, P2, and P3 for the same value of x5max. These equilibrium points can be stable or unstable, and the stability at the different steady states can be determined by calculating the eigenvalues of the matrix of the linearized model at points P1, P2, and P3. The value of x5max corresponding to the point M, (x5max)M, represents a global bifurcation. Thus, if the saturated flow is above this point the reactor behaviour Table 2 Dimensionless variables Dimensionless variable
Definition
Dimensionless variable
Definition
x1 x2 x3 x3r x4 x5 x5max x6 x20 x30 x40
Fos :t=Vs V=Vs Ca =Caor R:T=E R:Tr =E R:Tj =E Fj =Fjs Fjmax =Fjs F=Fos Cao =Caor R:To =E R:Tjo =E
x60 co c1* c11* c1 c2 c3 c4* c4 Kpd Kid Kd
Fo =Fos Vs :=Fos Vs..H.R. Caor Fos..cp.E c 1 =c 11 U:A=:cp :Fos Vs :Fjs =Fos :Vj :cp Vs :=j cpj Vj c4*.c2 Kp : To =Fjs Ki :ðFos =Vs Þ Kpd =x30
Subscript s denotes steady-state values; Caor is the initial value of reactant concentration in the reactor.
M. Pe´rez, P. Albertos / Journal of Process Control 14 (2004) 51–59
is similar to the unconstrained case, and the control system can drive the reactor to the set point. Also, if x5max is lower than (x5max)N, there is only one point of equilibrium at very high temperatures. Due to the low cooling water flow, the reactor temperature cannot reach an appropriated set point.
4. Self-oscillating conditions It is well known that a CSTR without control can present self-oscillations Vaganov [5,27]. From an industrial point of view, it is interesting to analyse the conditions under which the self-oscillations can arise in a controlled reactor. Let us study the appearance of self-oscillations by considering a constant reactor volume, without constrain in the cooling water flow, and applying P or PI control. In the following, the simplified dimensionless model of the reactor, 3, Eq. (9), will be used. The complete model 4, that is, considering the time constant of the jacket temperature, does not introduce qualitative differences and the computation burden increases. Let us first consider a P controller. In this case, the integral action is eliminated from Eqs. (9), taking the reset time equal to infinite. Considering that there is no disturbance in the inlet flow, 3 can be simplified putting x1=1, giving: dx2 ¼ x60 ðx20 x2 Þ c0 x2 e1=x3 d dx3 ¼ x60 ðx30 x3 Þ þ c1 x2 e1=x3 d c2 c3 x5 ðx3 x40 Þ c3 x5 þ c4 dx5 Kpd dx3 dx3 ¼ ¼ Kd d x30 d d
ð13Þ
Fig. 2. Dimensionless bifurcation plot when the cooling water is constrained. Parameters values: F0=1.7 m3/h, Ca0=8 kmol/m3, T0=Tj=21.7 C.
55
Integrating the last equation of (13), and substituting in the second one, the mathematical model of the Pcontrolled reactor is the following: dx2 1=x3 d ¼ x60 ðx20 x2 Þ c0 x2 e dx3 ¼ x60 ðx30 x3 Þ þ c1 x2 e1=x3 2 : d c2 c3 ð1 þ Kd ðx3 x3r ÞÞ ðx3 x40 Þ c3 ð1 þ Kd ðx3 x3r ÞÞ þ c4 ð14Þ The equilibrium reactor temperature, x3*, can be computed as a function of the x60, x30, x20, Kd and x3r: c1 x60 x20 x30 ¼ x 3
x60 e1=x3 þ c0 c2 c3 1 þ Kd x 3 x3r x 3 x40 þ ð15Þ c3 1 þ Kd x 3 x3r þ c4 From Eq. (15), assuming a constant value for x60, Kd and x3r, the equilibrium reactor temperature, x3*, can be expressed as a function, x 3 ¼ f3 ðx30 ; x20 Þ, of the inlet flow temperature and the inlet concentration. Some plots are shown in Fig. 3. For some given input conditions x 30 ; x 20 the number of equilibrium points depends on the number of inter
sections between the line x30 ¼ x30 and the function
x3 ¼ f3 x30 ; x20 . For example, for x30 ¼ 0:03 and x 20 ¼ 1:5, (see Fig. 3), there are three equilibrium points, P1, P2, P3, where P1 and P3 are stable and P2 is unstable. For this x 20 , a bifurcation point appears when the vertical line, x30=constant, is tangent to the curve df ðx ;x Þ x 3 ¼ f3 x30 ; x 20 . Thus, by solving 3 dx3030 20 ¼ 1, the necessary condition to have a bifurcation point for different equilibrium temperature, x3*, can be computed. If there is a real solution, x~ 30 , the pair x~ 30 ; x 20 represents a
Fig. 3. Inlet dimensionless temperature x30 vs. equilibrium point x3* at different values of x20. Ca0r=3.92 kmol A/m3; Kp=0.0066 m3/h C; Tr=60.6 C; x60=1.
56
M. Pe´rez, P. Albertos / Journal of Process Control 14 (2004) 51–59
possible bifurcation point for x3*. Parameterised on x3*, the possible bifurcation points are plotted in Fig. 4, splitting the input plane (x30, x20) into the inner subspace S1, with three equilibrium points, and the outer subspace, S2, with just one equilibrium point. Moreover, by the analysis of the local behaviour of the reactor around the equilibrium point defined by x3*, the existence of oscillatory modes can be determined. In this way, using Eqs. (14), the eigenvalues of the linearized model at the equilibrium point can be computed, and the pairs ðx 30 ; x 20 Þ leading to a couple of pure imaginary roots (the real part equal to zero) can be obtained. At these points, a local Hopf bifurcation appears. Parameterised on x3*, the points leading to such a pair of eigenvalues are also plotted in Fig. 4, splitting the input plane (x30, x20) into the inner subspace S3 and the outer subspace, S4. Any point in the intersection of
Fig. 4. Lobe with self-oscillating zone and curve with cusp point for the inlet stream dimensionless concentration and temperature.
Fig. 5. Eigenvalues of the model 3 [Eq. (9)] versus the volumetric flow rate inlet stream F0.
S2 and S3 leads to a self-oscillating behaviour characterized by a limit cycle. A similar approach can be followed in the case of a PI. The main difference is that the last equation in (13) cannot be integrated and, therefore, the model of the reactor, similar to (9), is a third order model. This introduces additional complexity in the computation with x60, x30, x20, Kpd, K and x3r as parameters, but following the same steps than before, the influence of the PI parameters in the steady state behaviour can be analysed.
5. Chaotic behaviour In Section 3 it was shown that a cooling water flow above the value (x5max)M, (see Fig. 2), avoids a global bifurcation and consequently, the control system can drive the reactor to the set point temperature. The change in the reactor behaviour under the influence of the input variables (x60, x20, x30, x40) is now analysed. It is possible, although really tedious, to calculate the variation of the eigenvalues of the linearized model around an equilibrium point when an input changes, whereas the others remain constant. In Fig. 5, the variation of the real part of eigenvalues vs. the volumetric flow F0 are plotted. It is clear that for high values of x60 ¼ F0 =F0s , there are eigenvalues with positive real part, and therefore the equilibrium point becomes unstable. Consider the state space model 3, with an equilibrium point such that the matrix of the linearized equation at this point has a real negative eigenvalue l and a pair of complex eigenvalues j. where > 0. In this situation, the equilibrium point has a onedimensional stable manifold and a two-dimensional unstable manifold. Imposing the condition jlj > , and supposing that there is a homoclinic orbit which tends to the equilibrium point, then the Silnikov theorem asserts that every neighbourhood of the homoclinic orbit contains an infinite number of unstable periodic orbits [8,21]. The existence of the equilibrium point and the conditions required to get these type of eigenvalues for the linearized model are not difficult to sort out. Nevertheless the hypothesis on the existence of a Silnikov homoclinic orbit is usually very difficult to establish. Following some hints from the literature, and after a trial-and-error simulation process, for x60 > 1 (see Fig. 5) and x5max(x5max)M (see Fig. 2), we found the presence of a homoclinic orbit to the equilibrium point. The adopted values for the rest of variables and parameters are: To=Tjo=21.7 C, Cao=8 kmol A/m3, Caor=3.92 kmol A/m3, where Caor is the initial value of reactant concentration in the reactor (see Table 2). The controller parameters are tuned to Kp=0.094 m3/h C,
M. Pe´rez, P. Albertos / Journal of Process Control 14 (2004) 51–59
57
Ki=5 h. The homoclinic orbit for the reactor simplified model is shown in Fig. 6. In accordance with Silkinov’s theorem, the reactor has a chaotic behaviour. In order to test the presence of a strange attractor, it is necessary to increase the value of x5max to introduce a perturbation in the vector field around the homoclinic orbit. Taking x5max=5 (Fjmax=7 m3/h), the results of the simulation are shown in Fig. 7. With the previous values, the sensitive dependence on initial conditions has been corroborated simulating Eqs. (8) and (9) with two initial conditions x01 and x02 very close in value. Taking x02 ¼ x01 þ 1 108 x01 , Fig. 8 shows that after a dimensionless time of 17 units, both solutions are completely different. It is interesting to note that the chaotic regime does not appear in the out stream flow which is manipulated by the control valve CV1, even if it is not in an initial
steady condition. So, the PI controller drives the reactor volume to the equilibrium point without the chaotic oscillations appearing on the other variables, as shown in Fig. 9. Thus, due to the decoupling already mentioned, it is possible to obtain a reactor behaviour, in which some variables are in steady state and others are in regime of chaotic oscillations. Note that the process nonlinearities, the limitation of the cooling water flow through the control valve CV2, as well as the control system are responsible for this interesting behaviour, which can appear in an industrial environment [26]. Note that this behaviour is completely different from that reported in the Pellegrini’s paper [15], where the route to chaos is due to a period doubling cascade. Moreover, in our analysis, the random behaviour appears taking into account the combined action of the PI controller and the control valve saturation. But it is
Fig. 6. Homoclinic Silnikov orbit type for the model 3 [Eq. (9). Eigenvalues at the equilibrium point for 3: l=0.566; j=0.398 1.925j, jlj > .
Fig. 8. Sensitive dependence of the dimensionless reactor concentration on two very close initial conditions. Parameter values: Kp=0.094 m3/h C, Ki=5 h, Fjmax=7 m3/h, F0=1.7 m3/h.
Fig. 7. Strange attractor for the model 3. Maximum volumetric flow rate of cooling water Fjmax=7 m3/h.
Fig. 9. Regular and chaotic oscillations in the reactor with PI control. Parameter values: Kp=0.094 m3/h C, Ki=5 h, Fjmax=7 m3/h, F0=1.7 m3/h.
58
M. Pe´rez, P. Albertos / Journal of Process Control 14 (2004) 51–59
also interesting to note that the PI controller actions above selected could be determined following the classical approach of Ziegler–Nichols. Even in this case, the chaotic behaviour may appear as far as the inlet temperature and concentration are both inside the lobe. Remark 3. In the simulation setting a fourth order Runge–Kutta or fifth order Runge–Kutta–Fehlberg methods have been used, with a step interval between 0.003 and 0.0005 units of dimensionless time. Nevertheless, many of the complicated orbits that arise after the breaking of the homoclinic orbit are unstable and, consequently they are not visible in the numerical computations. Remark 4. Similar results have been obtained using the full reactor model, 4.
6. Conclusions From the results of the research presented in this paper, it has been corroborated that a CSTR with an exothermic first order irreversible reaction, and PI control can be operated using a strategy that allows steady states, self-oscillation dynamics as well as to reach a more complex behaviour. The results for the CSTR show that no external forcing is necessary in order to obtain chaotic behaviour, since the process nonlinearities, the control system and the limitation of the cooling water flow are the origin of the chaotic oscillations. The variation of the composition and temperature of the reactor inlet stream shows that the steady state can be modified without any change in the controller parameters. It has been shown that with only proportional control, if the inlet flow dimensionless concentration x20 and temperature x30 belong to a small zone in the parameter space x20–x30, the reactor reaches a selfoscillating state. The existence and possible computation of bifurcation points is shown. This behaviour, which is independent of the proportional gain of the controller, cannot be obtained with a linear model. The presence of bifurcations and self-oscillations is also verified as a function of the set point temperature and the values of the proportional and integral constant of the PI controller, which manipulates the cooling water flow. The presence of a homoclinic orbit at the equilibrium point, often referred as Silnikov orbit, which implies chaotic dynamic, has been verified. The cooling water valve saturation is crucial to reach this behaviour. A system with a Silnikov orbit is a very interesting case, because the chaotic behaviour appears at and around the parameter value generating the homoclinic orbit. Thus, in this paper, a new case of a chaotic oscillation in a CSTR with a very simple control strategy is shown.
The previously analysed phenomena should be known by the control engineer in order to either avoid or properly use them, depending on the process type. Using the procedures shown in this paper, a CSTR with more sophisticated mathematical model and control system can be used to develop experimental and theoretical studies, to investigate the stability, self-oscillation and chaotic behaviour of this kind of processes.
References [1] F. Teymour, W.H. Ray, The dynamic behavior of continuous solution polymerization reactors-IV. Dynamic stability and bifurcation analysis of an experimental reactor, Chem. Eng. Sci. 44 (9) (1989) 1967–1982. [2] F. Teymour, Dynamics of semibatch polymerization reactors: I. Theoretical analysis, A.I.Ch.E. J. 43 (1) (1997) 145–156. [3] M. Soroush, Nonlinear state-observer design with applications to reactors, Chem. Eng. Sci. 52 (3) (1997) 387–404. [4] G. Kevrekidis, R. Aris, L.D. Schmidt, The stirred tank forced, Chem. Eng. Sci. 41 (6) (1986) 1549–1560. [5] D.A. Vaganov, N.G. Samoilenko, V.G. Abranov, Periodic regimes of continuous stirred tank reactors, Chem. Eng. Sci. 33 (1978) 1133–1140. [6] J.C. Mankin, J.L. Hudson, Oscillatory and chaotic behavior of a forced exothermic chemical reaction, Chem. Eng. Sci. 39 (12) (1984) 1807–1814. [7] M.L. Cartwright, J.E. Littlewood, On linear differential equations of the second order, London Mathematical Society 20 (1945) 180–189. [8] J. Guckenheimer, P. Holmes, Nonlinear Oscillations, Dynamical Systems, and Bifurcations of Vector Fields, Springer-Verlag, New York, 1983. [9] A.J. Lichtenberg, M.A. Lieberman, Regular and Chaotic Dynamics, 2nd Edition, Springer-Verlag, New York, 1992. [10] R. Seydel, Practical Bifurcation and Stability Analysis From Equilibrium to Chaos, 2nd Edition, Springer-Verlag, New York, 1994. [11] W.L. Luyben, Process Modeling, Simulation and Control for Chemical Engineers, 2nd Edition, McGraw-Hill, New York, 1990. [12] G. Stephanopoulos, Chemical Process Control. An Introduction to Theory and Practice, Prentice-Hall, Englewood Cliffs, 1991. [13] T. Marlin, Process Control, McGraw-Hill, New York, 1995. [14] B. Ogunnaike, W. Ray, Chemical Control Process, Oxford University, 1994. [15] L. Pellegrini, G. Biardi, Chaotic behaviour of a controlled CSTR, Computers Chem. Engng. 14 (1990) 1237–1247. [16] M.J. Kurtz, G. Yan Zhu, M.A. Henson, Constrained ouput feedback Control of a multivariable polymerization reactor, IEEE Transactions Control System and Technology 8 (2000) 87– 97. [17] J. Alvarez-Ramirez, R. Suarez, R. Femat, Robust stabilization of temperature in continuous-stirred tank reactors, Chem. Eng. Sci. 52 (14) (1997) 2223–2230. [18] J. Alvarez-Ramirez, R. Femat, J. Gonzalez-Trejo, Robust control of a class of uncertain first-order systems with least prior knowledge, Chem. Eng. Sci. 53 (15) (1998) 2701–2710. [19] R. Femat, Chaos in a class of reacting systems induced by robust asymptotic feedback, Physica D 136 (2000) 193–204. [20] M. Dolnik, A.S. Banks, I.R. Epstein, Oscillatory chemical reaction in a CSTR with feedback control of flow rate, Journal Phys. Chem. A 101 (1997) 5148–5154.
M. Pe´rez, P. Albertos / Journal of Process Control 14 (2004) 51–59 [21] L.P. Silnikov, A case of the existence of a demensurable set of periodic motions, Soviet Mathematical Dokl. 6 (1965) 163– 166. [22] D.R. Coughanowr, Process Systems Analysis and Control, 2nd Edition, McGraw-Hill, New York, 1991. [23] W.H. Ray, Advanced Process Control, McGraw-Hill, New York, 1981. [24] J.B. Planeaux, K.F. Jensen, Bifurcation phenomena in CSTR
59
dynamics: a system with extraneous thermal capacitance, Chem. Eng. Sci. 41 (6) (1986) 1497–1523. [25] R. Aris, Elementary Chemical Reactor Analysis, Prentice-Hall, Englewood Cliffs, NJ, 1969. [26] R. Aris, N.R. Amundson, An analysis of chemical reactor stability and control, Chem. Eng. Sci. 7 (1958) 121–130. [27] Z. Kubickova, M. Kubicek, M. Marek, Feed-batch operation of stirred reactors, Chem. Eng. Sci. 42 (2) (1987) 327–333.
Self-oscillating and chaotic behaviour of a PI-controlled CSTR with control valve saturation Manuel Pe´reza,1, Pedro Albertosb,* a
Department of ‘‘Fı´sica, Ingenierı´a de Sistemas y Teorı´a de la Sen˜al’’, Escuela Polite´cnica Superior, Alicante University, Campus de San Vicente, Alicante, Spain b Department of Systems Engineering and Control, E.T.S.I.I., Universidad Polite´cnica de Valencia, Camino de Vera s/n, Valencia, Spain Received 21 May 2002; received in revised form 16 April 2003; accepted 26 April 2003
Abstract The time response of a PI-controlled continuous stirred tank reactor (CSTR), where a first order irreversible reaction occurs, is analysed in this paper. The parameters of the PI controller and the saturation of the cooling water temperature flow are used as parameters to show that self-oscillations and chaotic dynamics can appear. The conditions under which self-oscillation occurs are investigated taking into account heat generation, heat removal and the variation of eigenvalues of the linearized system at the equilibrium points, respectively. Also, the appearance of a Silnikov type homoclinic orbit, giving chaotic dynamics, is realised. The existence of a new set of strange attractors with PI control has been verified and the sensitivity to the initial conditions is also analysed. # 2003 Elsevier Ltd. All rights reserved. Keywords: Control of continuous reactor; Steady state; Self-oscillation; Nonlinear and chaotic behaviour
1. Introduction The study of a continuous stirred tank reactor CSTR where a simple irreversible reaction A ! B occurs is basic for the analysis and control of chemical reactors. In the early years, this simple reaction scheme was considered as an ‘‘academic work’’ by many people, especially from the industry. However, the dynamics and the control of many reactors can be treated as processes in which the CSTR model is applicable. Processes involving catalytic reactors and hydrolysis reactions are, among others, good examples. The issues about the dynamic behaviour of CSTRs have been treated in various papers, [1–3], and books, and the oscillatory dynamics of a stirred tank reactor with a single reaction has been investigated using new numerical procedures in Kevrekidis [4]. In the past, various researchers [5,6] investigated a single irreversible * Corresponding author. Tel.: +34-96-387-9570; fax: +34-96-3879579. E-mail addresses: [email protected] (M. Pe´rez), [email protected] (P. Albertos). 1 Fax: +34-96-590-3682. 0959-1524/03/$ - see front matter # 2003 Elsevier Ltd. All rights reserved. doi:10.1016/S0959-1524(03)00032-5
exothermic reaction in which the coolant temperature was varied, following sine and square waves. They showed that with the variation of a single parameter, i.e. the amplitude of the variation in the coolant temperature, the reactor changes from quasi-periodic to periodic behaviour, and through a sequence of period doubling, it bifurcates to chaos. No control systems were used in any of these previous research studies. This is probably due to the fact that non-isothermal CSTR can show three steady states, not always being stable. The reactor evolves to one of its stable steady states [2] and without control, it is impossible to reach any arbitrary desired set point. Forced non-linear oscillators have been considered in various contexts. For example, the autonomous van der Pol oscillator, which is taken from electronics, is treated in Cartwright [7] showing the appearance of chaotic behaviour. Chaos can also occur via periodic perturbations in the Duffing equation [8]. Other examples can be found in Lichtenberg [9] and Seydel [10]. The simplest proportional controllers have been used to control the tank level and the temperature by manipulating the flow of the liquid leaving the tank and the cooling water [11–14]. In the paper of Pellegrini [15],
52
M. Pe´rez, P. Albertos / Journal of Process Control 14 (2004) 51–59
in order to get the chaotic behaviour of a PI-controlled CSTR, the outlet temperature is assumed as the manipulated variable without considering the cooling process dynamics. Moreover, more attention is paid to the non linear dynamics of the CSTR whereas the effect of the controller is not fully investigated. Also, in these conditions, it seems to be easier to reach the chaotic behaviour [6]. More recently, other more sophisticated control strategies for polymerization reactors based on feedback linearization and non-linear predictive control [16], or state observers [3] have been reported. In the work of Alvarez-Ramirez [17,18] the problem of the robust stabilisation of a CSTR has been studied and Femat [19], using robust asymptotic feedback, has corroborated the existence of the chaotic behaviour. In other papers, the self-oscillation in a CSTR using the flow rate as manipulated variable has been studied [20]. The PI control strategy for a CSTR provides a good example to analyse the difficulties encountered in controlling a process with several steady states. In this paper, an exothermic first order irreversible reaction in a CSTR is considered. The control system is very simple and common in the industry being composed by two PI controllers. The first controller manipulates the flow rate of the liquid leaving the reactor as a function of the volume in the tank error. It does not influence the qualitative response of the reactor and it can be assumed that it works perfectly, keeping constant the reactor volume. A second PI controller manipulates the flow rate of cooling water to the jacket as a function of the error in the reactor’s temperature. The paper shows that self-oscillation and chaotic behaviours are possible without any external periodic excitation. Just with proportional feedback control the reactor may reach self-oscillations independently of the values of the proportional gain. Another typical situation, like the saturation in the refrigerating control valve, is considered. This saturation may lead to a chaotic behaviour. Simulation studies are used for corroborating the existence of a homoclinic orbit of a special type, the so called Silnikov orbit [21], which contains an infinite number of unstable periodic orbits, that is, a chaotic behaviour.
2. CSTR model equations Let us assume an exothermic reaction A ! B in a CSTR. The reaction heat is removed by a cooling jacket surrounding the reactor. A first-order reaction in reactant A, negligible heat losses and perfectly mixing are assumed [11–13,22,23]. The jacket water is assumed to be perfectly mixed and the mass of the metal walls is considered negligible, so the thermal inertia of the metal is not considered. The reactor model without control includes the following equations:
dV ¼ Fo F dt
ð1Þ
d ðV C a Þ ¼ Fo Cao F Ca V Ca eE=RT dt
ð2Þ
dðVhÞ ¼ ðFo ho FhÞ HVCa eE=RT UA T Tj dt h ¼ cp T ; ho ¼ cp To ð3Þ
j Vj
dhj ¼ j Fj hjo hj þ UA T Tj dt hj ¼ cpj Tj ; hjo ¼ cpj Tjo
ð4Þ
obtained by a total mass balance (1), a component A mass balance (2), an energy balance in the reactor (3), and an energy balance in the jacket (4). A description of the variables and parameter values is given in Table 1. In addition to Eqs. (1)–(4), there are two PI controllers for the reactor outflow, F, and the cooling water flow, Fj, as shown in Fig. 1. The first one is acting on the inlet reactor flow to keep the volume constant. Note that the reactor volume, V(t), via manipulation of the outlet stream, is not affected by the rest of variables. Thus, in order to simplify the analysis, a perfect control of this volume is assumed, keeping it constant and equal to Vs. On the other hand, the volume does affect the temperature, the concentration and, thus, the refrigerating equations, and could be considered as a disturbance. The second one has the following equation: dFj dT 1 þ ðTr TÞ ¼ Kp ð5Þ dt Ki dt where Tr is the set point (operating point) of the reactor temperature and Kp, Ki are the controller proportional and integral action. The manipulated variable Fj is constrained: 0 4 Fj 4 Fjmax
ð6Þ
where Fjmax is the maximum value of the cooling water flow through the control valve CV2 (Fig. 1). The set of Eqs. (1)–(5) determine a good reactor model [13]. The variables Fo, To, Cao, and Tjo can be considered as load disturbances whose effect must be minimized by the control system. Remark 1. It will be shown that under nominal operating conditions, that is, V(t)=Vs and without any flow disturbance [Fo(t)=Fos], yielding VðtÞ=Fo ðtÞ constant, the reactor chaotic behaviour cannot be achieved despite other possible disturbances, such as To, Cao, and Tjo.
M. Pe´rez, P. Albertos / Journal of Process Control 14 (2004) 51–59
Remark 2. The dynamic related to the jacket temperature is much faster than that related to the reactor temperature. Thus, the jacket temperature time constant is negligible and a simplified model can be derived. With this assumption, from Eq. (4), it yields
Tj ¼
53
j cpj Fj Tjo þ UA T j cpj Fj þ UA
ð7Þ
To get the reduced model, Tj should be replace by this expression in Eq. (3).
Table 1 Variables, nominal operating conditions and parameter values Variable
Description
Value
V F Fo Ca h ho hj hjo T Tj Fj Fos Vs Cao Caor To Vj a E U A Tjo R H cp cpj r rj Tr
Reactor volume (m3) Volumetric flow rate for the outlet stream (m3/h) Volumetric flow rate for the inlet stream (m3/h) Reactant concentration for the outlet stream (kmol A/m3) Enthalpy for the outlet liquid (kJ/kg) Enthalpy for the inlet liquid (kJ/kg) Enthalpy for the outlet cooling water Enthalpy for the inlet cooling water Reactor temperature ( C) Jacket temperature ( C) Volumetric flow rate of cooling water (m3/h) Volumetric flow rate for the inlet stream (m3/h) Steady state reactor volume (m3) Reactant concentration inlet stream (kmol A/m3) Initial reactant concentration (kmol A/m3) Inlet stream temperature ( C) Jacket volume (m3) Preexponential factor from Arrhenius law (h1) Activation energy (kJ/kmol) Overall heat transfer in the jacket kJ/(h m2 C) Heat transfer area (m2) Inlet stream cooling water temperature ( C) Perfect-gas constant kJ/kmol K Enthalpy of reaction (kJ/kmol) Heat capacity inlet and out streams (kJ/kg K) Heat capacity of cooling water(kJ/kg K) Density of the inlet and out streams (kg/m3) Density of cooling water (kg/m3) Set point temperature ( C)
1.13 1.36 8 3.92 21.7 0.085 7.08 1010 69,815 3065 23.22 21.7 8.314 69,815 3.13 4.18 800 1000 60.6
Fig. 1. Perfectly mixed CSTR with two PI controllers.
54
M. Pe´rez, P. Albertos / Journal of Process Control 14 (2004) 51–59
2.1. CSTR dimensionless model To generalise the mathematical model of the reactor, the volume V being taken as a constant disturbance, Eqs. (2)–(5) can be expressed in a dimensionless way, leading to the following state-space model: dx2 x60 ¼ ðx20 x2 Þ c0 x2 e1=x3 d x 1 c2 dx3 x60 ¼ ðx30 x3 Þ þ c1 x2 e1=x3 ðx3 x4 Þ d x x 1 1 dx4 ¼ c x ð x x Þ þ c ð x x Þ 3 5 40 4 4 3 4 4 : d dx K x pd 5 60 ¼ ðx30 x3 Þ þ c1 x2 e1=x3 d x30 x1 c2 x3 x3r ðx3 x4 Þ þ x1 Kid
If the cooling water flow, Fjs, is not limited, providing as much refrigeration as required, the equilibrium points are obtained from (8), for x1e=1 and x3e ¼ x3r ;
x6e ¼ x60
ð10Þ
But in practice, as already noted in Eq. (6), this flow is limited and an equilibrium can be reached under saturated conditions. If this is the case, Fjs=Fjmax and the cooling water flow cannot decrease the reactor temperature, thus the equilibrium temperature reaches a value, Te, greater than the set point, (Te > Tr), and the corresponding reactant concentration will be smaller. The new equilibrium conditions are x5e ¼ x5max ;
x6e ¼ x6o
ð11Þ
Taking into account Eqs. (9) at the equilibrium points (11), the following relationship between the steady-state values can be obtained: ð8Þ
A full description of the dimensionless variables is given in Table 2. The fourth order model 4, (8), represents the reactor dynamics. If remark 2 is applicable, using Eq. (7) in (8), the reduced order model 3 will be dx2 x60 1=x3 d ¼ x1 ðx20 x2 Þ c0 x2 e dx3 x60 1=x3 d ¼ x1 ðx30 x3 Þ þ c1 x2 e c2 c3 x5 ðx3 x40 Þ 3 : x1 ðc3 x5 þ c4 Þ dx 5 ¼ KPd x60 ðx x Þ þ c x e1=x3 30 3 1 2 d x30 x1 c2 c3 x5 ðx3 x40 Þ þ ðx3 x3r Þ x1 ðc3 x5 þ c4 Þ Kid ð9Þ
3. Analysis of the effect of the control valve saturation It is well known that an exothermic CSTR without control has multiple steady states [24–27]. This is due to the fact that the curve of generated heat and the line of removed heat versus temperature can intersect at three points. A detailed behaviour of a CSTR without control can be found in Luyben [11], Stephanopoulos [12], Marlin [13] and Aris [25]. To run the reactor at any desired temperature T, the heat removal curve or the heat generation lines must be modified by changing some system parameters or adding a feedback control system. If a PI controller is used, Eqs. (2)–(4) should be considered but adding Eq. (5).
c2 c3 x5max x6e ðx30 x3e Þ x3e x40 c3 x5max þ c4 x6e x20 c1 ¼0 þ ðx3e x40 Þ ðx6e e1=x3e þ c0 Þ
ð12Þ
Eq. (12) shows the variation of the equilibrium dimensionless temperature as a function of the maximum value of dimensionless cooling water flow x5max. The plot of the implicit function f(x3e,x5max)=0, Eq. (12), is shown in Fig. 2. This plot shows the possibility of three steady states corresponding to the points P1, P2, and P3 for the same value of x5max. These equilibrium points can be stable or unstable, and the stability at the different steady states can be determined by calculating the eigenvalues of the matrix of the linearized model at points P1, P2, and P3. The value of x5max corresponding to the point M, (x5max)M, represents a global bifurcation. Thus, if the saturated flow is above this point the reactor behaviour Table 2 Dimensionless variables Dimensionless variable
Definition
Dimensionless variable
Definition
x1 x2 x3 x3r x4 x5 x5max x6 x20 x30 x40
Fos :t=Vs V=Vs Ca =Caor R:T=E R:Tr =E R:Tj =E Fj =Fjs Fjmax =Fjs F=Fos Cao =Caor R:To =E R:Tjo =E
x60 co c1* c11* c1 c2 c3 c4* c4 Kpd Kid Kd
Fo =Fos Vs :=Fos Vs..H.R. Caor Fos..cp.E c 1 =c 11 U:A=:cp :Fos Vs :Fjs =Fos :Vj :cp Vs :=j cpj Vj c4*.c2 Kp : To =Fjs Ki :ðFos =Vs Þ Kpd =x30
Subscript s denotes steady-state values; Caor is the initial value of reactant concentration in the reactor.
M. Pe´rez, P. Albertos / Journal of Process Control 14 (2004) 51–59
is similar to the unconstrained case, and the control system can drive the reactor to the set point. Also, if x5max is lower than (x5max)N, there is only one point of equilibrium at very high temperatures. Due to the low cooling water flow, the reactor temperature cannot reach an appropriated set point.
4. Self-oscillating conditions It is well known that a CSTR without control can present self-oscillations Vaganov [5,27]. From an industrial point of view, it is interesting to analyse the conditions under which the self-oscillations can arise in a controlled reactor. Let us study the appearance of self-oscillations by considering a constant reactor volume, without constrain in the cooling water flow, and applying P or PI control. In the following, the simplified dimensionless model of the reactor, 3, Eq. (9), will be used. The complete model 4, that is, considering the time constant of the jacket temperature, does not introduce qualitative differences and the computation burden increases. Let us first consider a P controller. In this case, the integral action is eliminated from Eqs. (9), taking the reset time equal to infinite. Considering that there is no disturbance in the inlet flow, 3 can be simplified putting x1=1, giving: dx2 ¼ x60 ðx20 x2 Þ c0 x2 e1=x3 d dx3 ¼ x60 ðx30 x3 Þ þ c1 x2 e1=x3 d c2 c3 x5 ðx3 x40 Þ c3 x5 þ c4 dx5 Kpd dx3 dx3 ¼ ¼ Kd d x30 d d
ð13Þ
Fig. 2. Dimensionless bifurcation plot when the cooling water is constrained. Parameters values: F0=1.7 m3/h, Ca0=8 kmol/m3, T0=Tj=21.7 C.
55
Integrating the last equation of (13), and substituting in the second one, the mathematical model of the Pcontrolled reactor is the following: dx2 1=x3 d ¼ x60 ðx20 x2 Þ c0 x2 e dx3 ¼ x60 ðx30 x3 Þ þ c1 x2 e1=x3 2 : d c2 c3 ð1 þ Kd ðx3 x3r ÞÞ ðx3 x40 Þ c3 ð1 þ Kd ðx3 x3r ÞÞ þ c4 ð14Þ The equilibrium reactor temperature, x3*, can be computed as a function of the x60, x30, x20, Kd and x3r: c1 x60 x20 x30 ¼ x 3
x60 e1=x3 þ c0 c2 c3 1 þ Kd x 3 x3r x 3 x40 þ ð15Þ c3 1 þ Kd x 3 x3r þ c4 From Eq. (15), assuming a constant value for x60, Kd and x3r, the equilibrium reactor temperature, x3*, can be expressed as a function, x 3 ¼ f3 ðx30 ; x20 Þ, of the inlet flow temperature and the inlet concentration. Some plots are shown in Fig. 3. For some given input conditions x 30 ; x 20 the number of equilibrium points depends on the number of inter
sections between the line x30 ¼ x30 and the function
x3 ¼ f3 x30 ; x20 . For example, for x30 ¼ 0:03 and x 20 ¼ 1:5, (see Fig. 3), there are three equilibrium points, P1, P2, P3, where P1 and P3 are stable and P2 is unstable. For this x 20 , a bifurcation point appears when the vertical line, x30=constant, is tangent to the curve df ðx ;x Þ x 3 ¼ f3 x30 ; x 20 . Thus, by solving 3 dx3030 20 ¼ 1, the necessary condition to have a bifurcation point for different equilibrium temperature, x3*, can be computed. If there is a real solution, x~ 30 , the pair x~ 30 ; x 20 represents a
Fig. 3. Inlet dimensionless temperature x30 vs. equilibrium point x3* at different values of x20. Ca0r=3.92 kmol A/m3; Kp=0.0066 m3/h C; Tr=60.6 C; x60=1.
56
M. Pe´rez, P. Albertos / Journal of Process Control 14 (2004) 51–59
possible bifurcation point for x3*. Parameterised on x3*, the possible bifurcation points are plotted in Fig. 4, splitting the input plane (x30, x20) into the inner subspace S1, with three equilibrium points, and the outer subspace, S2, with just one equilibrium point. Moreover, by the analysis of the local behaviour of the reactor around the equilibrium point defined by x3*, the existence of oscillatory modes can be determined. In this way, using Eqs. (14), the eigenvalues of the linearized model at the equilibrium point can be computed, and the pairs ðx 30 ; x 20 Þ leading to a couple of pure imaginary roots (the real part equal to zero) can be obtained. At these points, a local Hopf bifurcation appears. Parameterised on x3*, the points leading to such a pair of eigenvalues are also plotted in Fig. 4, splitting the input plane (x30, x20) into the inner subspace S3 and the outer subspace, S4. Any point in the intersection of
Fig. 4. Lobe with self-oscillating zone and curve with cusp point for the inlet stream dimensionless concentration and temperature.
Fig. 5. Eigenvalues of the model 3 [Eq. (9)] versus the volumetric flow rate inlet stream F0.
S2 and S3 leads to a self-oscillating behaviour characterized by a limit cycle. A similar approach can be followed in the case of a PI. The main difference is that the last equation in (13) cannot be integrated and, therefore, the model of the reactor, similar to (9), is a third order model. This introduces additional complexity in the computation with x60, x30, x20, Kpd, K and x3r as parameters, but following the same steps than before, the influence of the PI parameters in the steady state behaviour can be analysed.
5. Chaotic behaviour In Section 3 it was shown that a cooling water flow above the value (x5max)M, (see Fig. 2), avoids a global bifurcation and consequently, the control system can drive the reactor to the set point temperature. The change in the reactor behaviour under the influence of the input variables (x60, x20, x30, x40) is now analysed. It is possible, although really tedious, to calculate the variation of the eigenvalues of the linearized model around an equilibrium point when an input changes, whereas the others remain constant. In Fig. 5, the variation of the real part of eigenvalues vs. the volumetric flow F0 are plotted. It is clear that for high values of x60 ¼ F0 =F0s , there are eigenvalues with positive real part, and therefore the equilibrium point becomes unstable. Consider the state space model 3, with an equilibrium point such that the matrix of the linearized equation at this point has a real negative eigenvalue l and a pair of complex eigenvalues j. where > 0. In this situation, the equilibrium point has a onedimensional stable manifold and a two-dimensional unstable manifold. Imposing the condition jlj > , and supposing that there is a homoclinic orbit which tends to the equilibrium point, then the Silnikov theorem asserts that every neighbourhood of the homoclinic orbit contains an infinite number of unstable periodic orbits [8,21]. The existence of the equilibrium point and the conditions required to get these type of eigenvalues for the linearized model are not difficult to sort out. Nevertheless the hypothesis on the existence of a Silnikov homoclinic orbit is usually very difficult to establish. Following some hints from the literature, and after a trial-and-error simulation process, for x60 > 1 (see Fig. 5) and x5max(x5max)M (see Fig. 2), we found the presence of a homoclinic orbit to the equilibrium point. The adopted values for the rest of variables and parameters are: To=Tjo=21.7 C, Cao=8 kmol A/m3, Caor=3.92 kmol A/m3, where Caor is the initial value of reactant concentration in the reactor (see Table 2). The controller parameters are tuned to Kp=0.094 m3/h C,
M. Pe´rez, P. Albertos / Journal of Process Control 14 (2004) 51–59
57
Ki=5 h. The homoclinic orbit for the reactor simplified model is shown in Fig. 6. In accordance with Silkinov’s theorem, the reactor has a chaotic behaviour. In order to test the presence of a strange attractor, it is necessary to increase the value of x5max to introduce a perturbation in the vector field around the homoclinic orbit. Taking x5max=5 (Fjmax=7 m3/h), the results of the simulation are shown in Fig. 7. With the previous values, the sensitive dependence on initial conditions has been corroborated simulating Eqs. (8) and (9) with two initial conditions x01 and x02 very close in value. Taking x02 ¼ x01 þ 1 108 x01 , Fig. 8 shows that after a dimensionless time of 17 units, both solutions are completely different. It is interesting to note that the chaotic regime does not appear in the out stream flow which is manipulated by the control valve CV1, even if it is not in an initial
steady condition. So, the PI controller drives the reactor volume to the equilibrium point without the chaotic oscillations appearing on the other variables, as shown in Fig. 9. Thus, due to the decoupling already mentioned, it is possible to obtain a reactor behaviour, in which some variables are in steady state and others are in regime of chaotic oscillations. Note that the process nonlinearities, the limitation of the cooling water flow through the control valve CV2, as well as the control system are responsible for this interesting behaviour, which can appear in an industrial environment [26]. Note that this behaviour is completely different from that reported in the Pellegrini’s paper [15], where the route to chaos is due to a period doubling cascade. Moreover, in our analysis, the random behaviour appears taking into account the combined action of the PI controller and the control valve saturation. But it is
Fig. 6. Homoclinic Silnikov orbit type for the model 3 [Eq. (9). Eigenvalues at the equilibrium point for 3: l=0.566; j=0.398 1.925j, jlj > .
Fig. 8. Sensitive dependence of the dimensionless reactor concentration on two very close initial conditions. Parameter values: Kp=0.094 m3/h C, Ki=5 h, Fjmax=7 m3/h, F0=1.7 m3/h.
Fig. 7. Strange attractor for the model 3. Maximum volumetric flow rate of cooling water Fjmax=7 m3/h.
Fig. 9. Regular and chaotic oscillations in the reactor with PI control. Parameter values: Kp=0.094 m3/h C, Ki=5 h, Fjmax=7 m3/h, F0=1.7 m3/h.
58
M. Pe´rez, P. Albertos / Journal of Process Control 14 (2004) 51–59
also interesting to note that the PI controller actions above selected could be determined following the classical approach of Ziegler–Nichols. Even in this case, the chaotic behaviour may appear as far as the inlet temperature and concentration are both inside the lobe. Remark 3. In the simulation setting a fourth order Runge–Kutta or fifth order Runge–Kutta–Fehlberg methods have been used, with a step interval between 0.003 and 0.0005 units of dimensionless time. Nevertheless, many of the complicated orbits that arise after the breaking of the homoclinic orbit are unstable and, consequently they are not visible in the numerical computations. Remark 4. Similar results have been obtained using the full reactor model, 4.
6. Conclusions From the results of the research presented in this paper, it has been corroborated that a CSTR with an exothermic first order irreversible reaction, and PI control can be operated using a strategy that allows steady states, self-oscillation dynamics as well as to reach a more complex behaviour. The results for the CSTR show that no external forcing is necessary in order to obtain chaotic behaviour, since the process nonlinearities, the control system and the limitation of the cooling water flow are the origin of the chaotic oscillations. The variation of the composition and temperature of the reactor inlet stream shows that the steady state can be modified without any change in the controller parameters. It has been shown that with only proportional control, if the inlet flow dimensionless concentration x20 and temperature x30 belong to a small zone in the parameter space x20–x30, the reactor reaches a selfoscillating state. The existence and possible computation of bifurcation points is shown. This behaviour, which is independent of the proportional gain of the controller, cannot be obtained with a linear model. The presence of bifurcations and self-oscillations is also verified as a function of the set point temperature and the values of the proportional and integral constant of the PI controller, which manipulates the cooling water flow. The presence of a homoclinic orbit at the equilibrium point, often referred as Silnikov orbit, which implies chaotic dynamic, has been verified. The cooling water valve saturation is crucial to reach this behaviour. A system with a Silnikov orbit is a very interesting case, because the chaotic behaviour appears at and around the parameter value generating the homoclinic orbit. Thus, in this paper, a new case of a chaotic oscillation in a CSTR with a very simple control strategy is shown.
The previously analysed phenomena should be known by the control engineer in order to either avoid or properly use them, depending on the process type. Using the procedures shown in this paper, a CSTR with more sophisticated mathematical model and control system can be used to develop experimental and theoretical studies, to investigate the stability, self-oscillation and chaotic behaviour of this kind of processes.
References [1] F. Teymour, W.H. Ray, The dynamic behavior of continuous solution polymerization reactors-IV. Dynamic stability and bifurcation analysis of an experimental reactor, Chem. Eng. Sci. 44 (9) (1989) 1967–1982. [2] F. Teymour, Dynamics of semibatch polymerization reactors: I. Theoretical analysis, A.I.Ch.E. J. 43 (1) (1997) 145–156. [3] M. Soroush, Nonlinear state-observer design with applications to reactors, Chem. Eng. Sci. 52 (3) (1997) 387–404. [4] G. Kevrekidis, R. Aris, L.D. Schmidt, The stirred tank forced, Chem. Eng. Sci. 41 (6) (1986) 1549–1560. [5] D.A. Vaganov, N.G. Samoilenko, V.G. Abranov, Periodic regimes of continuous stirred tank reactors, Chem. Eng. Sci. 33 (1978) 1133–1140. [6] J.C. Mankin, J.L. Hudson, Oscillatory and chaotic behavior of a forced exothermic chemical reaction, Chem. Eng. Sci. 39 (12) (1984) 1807–1814. [7] M.L. Cartwright, J.E. Littlewood, On linear differential equations of the second order, London Mathematical Society 20 (1945) 180–189. [8] J. Guckenheimer, P. Holmes, Nonlinear Oscillations, Dynamical Systems, and Bifurcations of Vector Fields, Springer-Verlag, New York, 1983. [9] A.J. Lichtenberg, M.A. Lieberman, Regular and Chaotic Dynamics, 2nd Edition, Springer-Verlag, New York, 1992. [10] R. Seydel, Practical Bifurcation and Stability Analysis From Equilibrium to Chaos, 2nd Edition, Springer-Verlag, New York, 1994. [11] W.L. Luyben, Process Modeling, Simulation and Control for Chemical Engineers, 2nd Edition, McGraw-Hill, New York, 1990. [12] G. Stephanopoulos, Chemical Process Control. An Introduction to Theory and Practice, Prentice-Hall, Englewood Cliffs, 1991. [13] T. Marlin, Process Control, McGraw-Hill, New York, 1995. [14] B. Ogunnaike, W. Ray, Chemical Control Process, Oxford University, 1994. [15] L. Pellegrini, G. Biardi, Chaotic behaviour of a controlled CSTR, Computers Chem. Engng. 14 (1990) 1237–1247. [16] M.J. Kurtz, G. Yan Zhu, M.A. Henson, Constrained ouput feedback Control of a multivariable polymerization reactor, IEEE Transactions Control System and Technology 8 (2000) 87– 97. [17] J. Alvarez-Ramirez, R. Suarez, R. Femat, Robust stabilization of temperature in continuous-stirred tank reactors, Chem. Eng. Sci. 52 (14) (1997) 2223–2230. [18] J. Alvarez-Ramirez, R. Femat, J. Gonzalez-Trejo, Robust control of a class of uncertain first-order systems with least prior knowledge, Chem. Eng. Sci. 53 (15) (1998) 2701–2710. [19] R. Femat, Chaos in a class of reacting systems induced by robust asymptotic feedback, Physica D 136 (2000) 193–204. [20] M. Dolnik, A.S. Banks, I.R. Epstein, Oscillatory chemical reaction in a CSTR with feedback control of flow rate, Journal Phys. Chem. A 101 (1997) 5148–5154.
M. Pe´rez, P. Albertos / Journal of Process Control 14 (2004) 51–59 [21] L.P. Silnikov, A case of the existence of a demensurable set of periodic motions, Soviet Mathematical Dokl. 6 (1965) 163– 166. [22] D.R. Coughanowr, Process Systems Analysis and Control, 2nd Edition, McGraw-Hill, New York, 1991. [23] W.H. Ray, Advanced Process Control, McGraw-Hill, New York, 1981. [24] J.B. Planeaux, K.F. Jensen, Bifurcation phenomena in CSTR
59
dynamics: a system with extraneous thermal capacitance, Chem. Eng. Sci. 41 (6) (1986) 1497–1523. [25] R. Aris, Elementary Chemical Reactor Analysis, Prentice-Hall, Englewood Cliffs, NJ, 1969. [26] R. Aris, N.R. Amundson, An analysis of chemical reactor stability and control, Chem. Eng. Sci. 7 (1958) 121–130. [27] Z. Kubickova, M. Kubicek, M. Marek, Feed-batch operation of stirred reactors, Chem. Eng. Sci. 42 (2) (1987) 327–333.