Method of determination of steady-state diagrams of chemical reactors

Chemical Engineering Science 55 (2000) 4291}4295

Method of determination of steady-state diagrams of chemical reactors Marek Berezowski 
* Cracow University of Technology, Institute of Chemical Engineering, 31-155 Krako& w, ul. Warszawska 24, Poland
Polish Academy of Sciences, Institute of Chemical Engineering, 44-100 Gliwice, ul. Baltycka 5, Poland Received 1 September 1999; received in revised form 7 December 1999; accepted 22 February 2000

Abstract The determination of pro"les and bifurcation diagrams of tubular reactors and CSTR cascade is described. This method makes use of the "ctitious dynamics and the continuation numerical process. It can be applied for di!erent reactor con"gurations and chemical reactions. In this method the optimization procedure to calculate the boundary conditions is not necessary.  2000 Elsevier Science Ltd. All rights reserved. Keywords: Bifurcation diagram; Steady state; Chemical reactor

1. Introduction One of the main problems in the process of design and theoretical analysis of chemical plants is the problem of steady states of a given system. In the case when a plant comprises elements described by non-linear di!erential equations with boundary conditions, the solution of steady-states problem (numerically or analytically) becomes di$cult and requires often the application of sophisticated convergence procedures. The problem becomes particularly troublesome when a system is characterized by the so-called multiplicity of steady states of a complex structure. This phenomenon occurs, "rst of all, in various types of chemical reactors and in systems comprising the latter. Quite often, not only the determination of a single steady state of a given system, but a full range of the occurrence of steady states is required, i.e. of the so-called bifurcation diagram, which all the more complicates the numerical calculations. The present paper presents a relatively simple, general method of the determination of bifurcation diagrams of steady-state systems comprising chemical reactors. The method is of general character and does not require either discretization of the model or application of iter* Corresponding author. Polish Academy of Sciences, Institute of Chemical Engineering, 44-100 Gliwice, ul. Baltycka 5, Poland. Fax: #0048-32-231-0318. E-mail address: [email protected] (M. Berezowski).

ative adjustment of boundary conditions. The method enables, in one computational course, to determine the complete bifurcation diagram of steady states. The application of this method has been illustrated by three examples. The "rst/second of these examples concerns a non-adiabatic tubular heterogeneous/homogeneous chemical reactor with longitudinal dispersion of heat and mass. Moreover, the reactor may have an additional external heat and/or mass coupling accomplished by recycle and/or heat exchanger. In the reactor an arbitrary number of r independent chemical reactions of arbitrary kinetics may take place (Berezowski, Ptaszek & Jacobson, 2000). The third example concerns a cascade of non-adiabatic CSTR. The cascade may be composed of an arbitrary number of CSTR and may have an external heat and/or mass coupling accomplished by recycle and/or external heat exchanger. Similarly as above, in the individual reactor an arbitrary kinetics may occur (Berezowski, Ptaszek, Grzywacz, Z ukowski, 1999). 2. General description of the method of determination of the bifurcation diagram of steady states of boundary value problem Let us assume that one has to do with a steady-state problem including, a.o., the boundary conditions of the form F [u (¸), u (1!¸)]"0, j"1, m; ¸"0 or 1, H

0009-2509/00/$ - see front matter  2000 Elsevier Science Ltd. All rights reserved. PII: S 0 0 0 9 - 2 5 0 9 ( 0 0 ) 0 0 0 5 8 - 0

(1)

4292

M. Berezowski / Chemical Engineering Science 55 (2000) 4291}4295

where u is the n-element vector of independent variables and m)n. Treating the variables u (¸) as variables of some "ctiH tious dynamic problem, we may introduce into Eq. (1) a "ctitious dynamics of the form l H

du (¸, q) H "F [u (¸, q), u (1!¸), q], j"1, m. H dq

(2)

System (2) constitutes an initial-value problem with condition u (¸, 0)"u (¸) and its steady-state solutions ful"l  identically the relationships with Eq. (1). The variable u (1!¸, q), appearing in Eq. (2), may be determined from remaining equations of the physical model (as a rule, the balance equations of the described processes). For that matter, these can be either di!erential equations (a continuous process): du G " (u ), i"1, n, ¸)1)1!¸ G d1

(3)

or algebraic ones (discrete process): u (¸!1)"t [u (¸)], i"1, n. (4) G G Also mixed relationships of Eqs. (3) and (4) type come into question. The parameters l should be selected in H such a way that system (2) be stable. As a rule, it is su$cient to assume l "$1. Thus, the method enables G one to determine both stable and unstable steady states of a real object. In order to determine the complete bifurcation diagram the method of continuation of the solutions of RHS of system (2) should be applied. As starting values in this method the solution obtained by the simulation as above should be assumed. As may be seen, the computational procedure is deprived of any iterations (Berezowski, 2000). 3. Example of the application of the method to the determination of diagram of steady states of a tubular reactor The general model of non-adiabatic tubular heterogeneous reactor with longitudinal dispersion of heat and mass and additional external coupling, in the steady state, can be presented as follows: Gas phase: da 1 da G" G # , i"1, r, G Pe dm dm + dH 1 dH $" $ #St(H !H ) I $ dm Pe dm &$ Solid phase:

The corresponding boundary conditions are: for m"0: 1 da G "a !f a (1), i"1, r, G + G Pe dm + 1 dH $ "H !f H (1), $ & $ Pe dm &$ dH I "0, dm

(8) (9) (10)

for m"1: da G "0, dm

dH $ "0, dm

dH I "0, i"1, r. dm

(11)

After substitutions: da u "a , u " G , i"1, r G G G dm

(12)

dH u "H , u " $ ,  $  dm

(13)

dH u " I ,  dm

(14)

u "H ,  I

model (5)}(11) may be transformed as follows: du G "u , i"1, r, G dm

(15)

du G "Pe (u ! ), i"1, r + G G dm

(16)

du  "u ,  dm

(17)

du  "Pe [u !St(u !u )], &$    dm

(18)

du  "u ,  dm

(19)





du GP  "Pe St(u !u )#d(u !H )! . &I    & G dm G The boundary conditions take the form:

(20)

for m"0: (5) (6)

1 dH GP I # !St(H !H )!d(H !H ). 0" G I $ I & Pe dm &I G (7)

1 u "u !f u (1), i"1, r, G + G Pe G + 1 u "u !f u (1),   &  Pe &$ u "0,  for m"1: u "0, u "0, u "0, i"1, r. G  

(21) (22) (23)

(24)

M. Berezowski / Chemical Engineering Science 55 (2000) 4291}4295

4293

Let us consider the "ctitious dynamics concerning the state variables: u (1), u (1) and u (1), connected with the G   boundary conditions (21)}(23) in the form of following di!erential equations: 1 du (1, q) G " u (0, q)!u (0,q)#f u (1, q), G G + G dq Pe + i"1, r, (25)

l G

du (1, q) 1  " u (0, q)!u (0, q)#f u (1, q),   &  dq Pe &$ du (1, q)  l "u (0, q).   dq l 

(26) (27)

It may clearly be seen that the steady state of the "ctitious model as above ful"ls identically the steady-state solutions of reactor's model. In order to perform the computer simulation of the solution of system (25)}(27) at an arbitrary moment q, this system should be integrated. The variables: u (0), u (0), u (0), u (0), u (0), u (0) are G G     determined by backward integration of Eqs. (11)}(16) with the initial conditions: u (1), u (1), u (1) and Eq. (24). G   Carrying out in the mentioned way the computations, one obtains, after a certain simulation time, a single steady-state solution (steady-state temperature * and conversion pro"les) of reactor's model Eqs. (5)}(11). In order to obtain the complete bifurcation diagram, as a function of a selected bifurcation parameter, one should successively apply the method of continuation (Doedel, et al., 1997) of solution of RHS of system (25)}(27). As starting values in this method one should assume the solution obtained by the simulation as above. Thus, the whole computational procedure is deprived of any iterations and leads to a quick determination of the diagram. The method described is of general character in such a manner that the introduction of additional coupling, in the form of additional boundary conditions, does not a!ect the way of obtaining the solution. In order to illustrate the above method two bifurcation diagrams of a tubular reactor will be determined. The "rst one deals with the non-adiabatic heterogeneous tubular reactor with longitudinal dispersion of heat and mass in the #uid phase and heat conduction in the solid phase. Furthermore, the reactor's recycle should be allowed for. A single chemical reaction APB type of arbitrary order:




"(1!f )Da(1!a)L exp c



bH I , 1#bH I

(where f"f "f ) (28) + & has been assumed. Assuming the following values of the parameters: Da"0.15, H "!0.06, c"15, n"1.5, & b"2, Pe "100, Pe "100, Pe "100, f"0.7, + &$ &I d"0.9 the diagram has been obtained (Fig. 1). The physical sense concerns, evidently, the interval St'0,

Fig. 1. Bifurcation diagram of steady states of tubular heterogeneous reactor.

Fig. 2. Bifurcation diagram of steady states of tubular homogeneous reactor.

but from the formal view point the whole loop has been shown. The second diagram concerns the homogeneous chemical reactor with longitudinal dispersion of heat and mass, as well as recycle. Assuming: Da"0.15, c"15, n"1.5, b"2, Pe "300, Pe "300, f"0.9, d"0.3 + & and St"R the diagram of hysteresis type has been constructed (Fig. 2).

4. Example of application of the method to the determination of diagrams of steady states of a CSTR cascade Let us assume that we have to do with a cascade composed of N non-adiabatic tank reactors (CSTR), in which r independent chemical reactions take place. The cascade is coupled through recycle and/or external heat exchanger. The model of such a system can be described, in the steady state, by the following recurrence sequence: a "a !

, ,\H G ,>\H G ,>\H G JP H "H !

,\H ,>\H ,>\HJ J !d [H !H ], ,>\H & ,>\H

(29)

(30)

4294

M. Berezowski / Chemical Engineering Science 55 (2000) 4291}4295

where i"1, r; j"1, N!1. The boundary conditions resulting from the couplings: a "f a , H "f H , i"1, r. G + ,G G & ,

(31)

Similarly as in the case of a tubular reactor let us introduce into considerations a "ctitious dynamics concerning state variables: a , H , (i"1, r). Substituting Eq. (31) ,G ,G into Eqs. (29) and (30) one obtains a "ctitious dynamic model in the form: l +G l &

da (q) ,G "f a (q)!a (q)# , i"1, r + ,G G G dq

(32)

dH (q) , "f H (q)!H (q) & ,  dq

Notation

JP # #d [H !H (q)]. J  &  J

(33)

It is clearly seen that the steady state of the above dynamic model ful"ls identically the steady-state solution of the model of CSTR cascade Eqs. (29)}(31). Hence, by integration of systems (32) and (33) one obtains, after a certain simulation time, a solution being arbitrarily close to the steady-state one. This solution, in turn, constitutes a starting point for computations by the continuation method (Doedel et al., 1997) aimed at the determination of the complete bifurcation diagram of steady states of the cascade. Strictly speaking, one continuation course enables one to determine the diagrams of all N reactors composing the cascade. On the other hand, the number of equations requiring integration is independent of cascade's length. Similarly as in the previous case, the method under consideration is of general character due to the form of boundary conditions as well as the number and type of chemical reactions. For the sake of illustration of the above method, the steady-state diagrams of reactors composing a cascade of N"100 CSTR, have been determined. The cascade is coupled by means of recycle and in individual reactors the reactions of APB occur with the kinetics:






bH H

"(1!f )Da(1!a )L exp c , H H 1#bH H j"1, N; (where f"f "f ). + &

Fig. 3. Bifurcation diagrams of steady states selected CSTR from cascade.

Da f n N Pe & Pe &$ Pe &I Pe + r St

DamkoK hler number feedback coe$cient order of reaction number of reactors in cascade Peclet number concerning heat dispersion in homogeneous reactor "Pe Pe /(Pe #Pe ) &I &$ &I &$ Peclet number concerning heat dispersion in #uid phase Peclet number concerning heat conduction in solid phase Peclet number concerning mass dispersion number of chemical reactions Stanton number

Greek letters a b c d H l m q

conversion degree dimensionless number related to adiabatic temperature rise dimensionless activation energy dimensionless heat transfer coe$cient dimensionless temperature real number dimensionless position coordinate along reactor dimensionless time dimensionless kinetic function

Subscripts (34)

Assuming the following values of the parameters: f"0.1, d"0.05, H "0, N"100, c"15, b"1, n"1.5, the & diagrams of the 1st, 5th, 10th, 50th and 100th reactor are presented (Fig. 3). It is notewothy that a number of stationary states, available in the ith CSTR, equals 3G (provided that a number of such states in a reactor outside a cascade equals 3).

F H k M

#uid phase heat, heat exchanger solid phase mass

References Berezowski, M. (2000). http://www.iich.gliwice.pl/users/mberez/www/ diagram.

M. Berezowski / Chemical Engineering Science 55 (2000) 4291}4295 Berezowski, M., Ptaszek, P., Grzywacz, R., & Z ukowski, W. (1999). Theoretical analysis of static and dynamic behaviour in systems based on a cascade of reactors. Inzynieria Chemiczna i Procesowa, 20, 185}207 (in Polish). Berezowski, M., Ptaszek, P., & Jacobsen, E. W. (2000). Dynamics of heat-integrated pseudohomogeneous tubular reactors with axial dispersion. Chemical Engineering and Processing, 39, 181}188.

4295

Doedel, E. J., Champneys, A. R., Fairgrieve, T. F, Kuznetsov, Y. A., Sandstede, B., & Wang, X. (1997). AUTO97: Continuation and bifurcation software for ordinary di!erential equations. Technical Report, Computational Mathematics Laboratory, Concordia University. Also http://indy.cs.concordia.ca/auto.