Combustion, Explosion, and Shock Waves, Vol. 37, No. 6, pp. 634 640, 2001
Parametric Analysis of the Models
of a Stirred Tank Reactor and a Tube Reactor
V. I. Bykov1 and S. B. Tsybenova2 UDC 536.46
Translated from Fizika Goreniya i Vzryva, Vol. 37, No. 6, pp. 22 29, November December, 2001.
Original article submitted August 15, 2000.
A parametric analysis is performed for a basic model of the combustion theory with
dimensional parameters. The model describes the dynamics of a first-order exother-
mic reaction in a well-stirred tank reactor. Parametric dependences of steady states
on dimensional parameters, curves of multiplicity and neutrality of steady states, and
parametric and phase portraits of the system are constructed. Regions of multiple
steady states and self-excited oscillations and a region of technological safety are dis-
tinguished. Combustion regimes in a continuous stirred tank reactor are compared
to those in a tube reactor.
INTRODUCTION PARAMETRIC ANALYSIS
OF A DIMENSIONAL MODEL
The authors [1, 2] performed a parametric analy-
The mathematical model of a spherical continuous
sis for models with dimensionless parameters and ob-
stirred tank reactor (Frank-Kamenetskii reactor) has
tained results of methodical and theoretical interest.
the form
However, in modeling specific processes and reactors,
dX
it is necessary to analyze models with parameters pos-
V = -V k(T )X + q(X0 - X) = f(X, T ),
dt
sessing a definite physicochemical meaning by numerical
and qualitative methods. Therefore, a parametric anal-
dT
ysis should also be performed for specific mathemati- cpÁV = (-"H)V k(T )X (1)
dt
cal models that correspond to the geometrical, thermal,
and physical characteristics of a reactor considered.
+ qcpÁ(T0 - T ) + hS(Tw - T ) = g(X, T ),
In the present paper, we illustrate the parametric
where V [cm3] is the volume of the reactor, k(T ) =
analysis procedure using a dimensional model of a con-
k0 exp(-E/RT ) [sec-1] is the reaction rate constant,
tinuous stirred tank reactor, where only one first-order
k0 is the preexponent, X and X0 [mole/cm3] are the
exothermic reaction occurs. The relationship between
current and inlet concentrations of the reagent, respec-
the results of the parametric analysis of dimensionless
tively, T and T0 [K] are the current and inlet tem-
and dimensional models is established. Moreover, we
peratures of the mixture, q [cm3/sec] is the volume
compare the dynamic properties of a continuous stirred
flow rate, Á [mole/cm3] is the density of the mixture,
(spherical) tank reactor and continuous tube (cylindri-
cp [J/(mole · K)] is the specific heat of the reactive mix-
cal) reactor characterized by similar physical and tech-
ture, S [cm2] is the area of the heat-exchange surface,
nological parameters in which the same chemical reac-
h [J/(cm2· sec · K)] is the coefficient of heat transfer
tion occurs.
through the reactor wall, (-"H) [J/mole] is the ther-
mal effect of the reaction, t [sec] is the astronomical
time, and Tw [K] is the temperature of the reactor walls.
1
Institute of Computational Modeling,
A parametric analysis of the dimensional model (1)
Siberian Division, Russian Academy of Sciences,
is performed according to the procedure applied to di-
Krasnoyarsk 660036.
2 mensionless models [1, 2].
Krasnoyarsk State Technical University,
Krasnoyarsk 660074.
634 0010-5082/01/3706-0634 $25.00 © 2001 Plenum Publishing Corporation
Parametric Analysis of the Models of a Stirred Tank Reactor and a Tube Reactor 635
Steady States. Equating the right sides of sys-
tem (1), we obtain the stationary conditions
-V k(T )X + q(X0 - X) = 0,
(2)
(-"H)V k(T )X + qcpÁ(T0 - T ) + hS(Tw - T ) = 0,
which imply that
qcpÁ(T0 - T ) + hS(Tw - T )
X = + X0. (3)
q(-"H)
Inserting (3) into the second equation in (2), we obtain
the equation for determining the steady-state tempera-
tures:
V k(T )
(-"H)V k(T )X0 + qcpÁ(T0 - T )
q
+hS(Tw - T ) + qcpÁ(T0 - T )
+ hS(Tw - T ) = 0. (4)
As the basic parameters of model (1), we use quan-
tities corresponding to the experimental data on oxida-
tion of hydrocarbons in a continuous stirred tank [3, 4].
In the parametric analysis, some parameters given in
Table 1 were slightly varied. The varied parameters in-
clude the inlet temperature and concentration and the
heat-transfer coefficient.
Parametric Dependences. Following the gen-
eral procedure of parametric analysis [2], from (4) we Fig. 1. Steady-state temperature versus the param-
eters T0 and X0 for model (1).
obtain the dependences inverse to the desired depen-
dences of the steady-state temperature and concentra-
tion on the parameters of the system:
variation in various parameters of the reactor consid-
1
ered. Parametric dependences of the steady-state con-
X0(T ) = cpÁ(T - T0)(q + V k(T ))
(-"H)V k(T )
centrations X are constructed with the use of relation
(3) where T is replaced by the corresponding depen-
V k(T )
dences T (X0), T (T0), etc.
+ hS(T - Tw) 1 + , (5)
q
For multiple steady states, the parametric depen-
dence of the steady-state temperature (for example,
1
T0(T ) = T - (-"H)V k(T )X0 curve 2 in Fig. 1a) has three branches of low, high, and
cpÁ(q + V k(T ))
intermediate temperatures. Transition from one branch
to another is performed by jumps. Upon attainment
V k(T )
- hS(T - Tw) 1 + . (6)
of a certain critical state, the temperature of the reac-
q
tor changes abruptly, i.e., its  ignition or  extinction
Similarly to (5) and (6), one can easily obtain de-
occurs.
pendences on other parameters: V (T ), S(T ), Tw(T ),
Stability of Steady States. Stability of steady
h(T ), Á(T ), etc. Figure 1 shows some parametric depen-
states is determined by the roots of the characteristic
dences T (X0) and T (T0) obtained by formulas (5) and
equation
(6). These curves are characterized by a region of multi-
2 - Ã + " = 0, (7)
valued solutions, which corresponds to multiple steady
states. With a second parameter varied (curves 3 5 in where à = a11+a22, " = a11a22-a12a21, and aij are the
Fig. 1b), the region of multiple steady states changes elements of the Jacobi matrix of the right sides of sys-
significantly or even disappears (curve 1). These para- tem (1), calculated in the steady states: a11 = "f/"X,
metric dependences allow one to estimate the conditions a12 = "f/"T , a21 = "g/"X, and a22 = "g/"T . The
under which hysteresis and critical phenomena exist for steady state is unstable if " < 0. If " > 0, the steady
636 Bykov and Tsybenova
TABLE 1
Basic Parameters of Model (1)
Physical parameters Values
Inlet temperature, T0 583 K
Inlet concentration of the substance, X0 1.4 · 10-6 mole/cm3
Density of the mixture, Á 3 · 10-5 mole/cm3
Specific heat of the mixture, cp 2166.2 J/(mole · K)
Volume of reactor, V 588.75 cm3
Area of the reactor surface, S 471 cm2
Coefficient of heat transfer, h 8 · 10-3 J/(cm2 · sec · K)
Reaction heat, (-"H) 4.87 · 106 J/mole
Volume flow rate, q 180 cm3/sec
Temperature of the reactor walls, Tw 583 K
Activation energy, E 1.89 · 105 J/mole
Universal gas constant, R 8.31 J/(mole · K)
Preexponent, k0 2.3 · 1015 sec-1
Å„Å‚
X0 = X0(T, T0),
state is stable for à < 0 and unstable for à > 0. Con-
ôÅ‚
ôÅ‚
ôÅ‚
A
òÅ‚
sequently, for " = 0 or à = 0, the type of stability of
T0 = T - (T - Tw)
L"(X0, T0): (8)
q
the steady states and their number change, i.e., a bifur-
ôÅ‚
ôÅ‚
(q + V k(T ))(q + A)
ôÅ‚
cation occurs. Therefore, the conditions that " and Ã
ół
+ ,
2
Eq2/RT
vanish yield important information on the number and
type of steady states. Combined with the stationary
L"(T0, Á):
condition (4), these equalities allow one to construct bi-
furcation curves for various combinations of parameters:
Å„Å‚
1
ôÅ‚
curves of multiplicity (" = 0) and neutrality (Ã = 0) of
ôÅ‚
T0(T, Á) = T - (-"H)
ôÅ‚
ôÅ‚
cpÁ(q + V k(T ))
ôÅ‚
steady states.
ôÅ‚
ôÅ‚
òÅ‚
Bifurcation Curves. For model (1), as in the
V k(T )
× V k(T )X0 - hS(T - Tw) 1 + ,
case of dimensionless models [1], we write equations for
ôÅ‚
ôÅ‚
q
ôÅ‚
ôÅ‚
constructing local bifurcation curves L" (" = 0) and
ôÅ‚
E qk(T )(-"H)X0 hS
ôÅ‚
ôÅ‚
ół Á(T ) = - ,
LÃ (Ã = 0) in various planes. Combining, for example,
2
RT cp(q + V k(T ))(k(T ) + q/V ) qcp
the equation "(X0, V ) = 0 with (4), one can obtain the
where A = hS/cpÁ.
dependence X0(T, V ), which, together with (5), defines
Neutrality Curves:
a multiplicity curve in the parameter plane (X0, V ).
Multiplicity Curves:
LÃ(X0, T0):
Å„Å‚
ôÅ‚ X0 = X0(T, T0),
L"(X0, V ):
ôÅ‚
ôÅ‚
ôÅ‚
A 1
òÅ‚
T0(T ) = T - - (T - Tw)
q B
ôÅ‚
ôÅ‚
1
ôÅ‚
ôÅ‚
ół - (2q + V k(T )),
Å„Å‚
qB
1
ôÅ‚
ôÅ‚
X0(T, V ) = cpÁ(T - T0)
ôÅ‚
ôÅ‚
(9)
(-"H)V k(T )
ôÅ‚
ôÅ‚
ôÅ‚
ôÅ‚
òÅ‚
Å„Å‚
V k(T )
X0 = X0(T, V ),
ôÅ‚
× (q + V k(T )) + hS(T - Tw) 1 + ,
ôÅ‚
òÅ‚
ôÅ‚ 1
q
ôÅ‚
ôÅ‚ V = (A(B(T - Tw) - 1)
ôÅ‚
LÃ(X0, V ):
ôÅ‚
k(T )
ôÅ‚
ôÅ‚
E q2(T - T0) + qA(T - Tw) q
ôÅ‚
ôÅ‚
ôÅ‚
ół
ół V (T ) = - ,
2
+ q(B(T - T0) - 2)) ,
RT k(T )(q + A) k(T )
Parametric Analysis of the Models of a Stirred Tank Reactor and a Tube Reactor 637
Fig. 3. Parametric portrait for model (1) in the
(X0, V ) plane for h = 0.04 J/(cm2 · sec · K) and
T0 = 583 K.
Fig. 2. Bifurcation curves of multiplicity L" (a) and
neutrality LÃ (b) for model (1) in the (X0, V ) plane.
LÃ(T0, Á):
Fig. 4. Phase portrait of system (1) for the unique
Å„Å‚
and unstable steady state for V = 1500 cm3, h =
ôÅ‚ T0 = T0(T, Á),
ôÅ‚
ôÅ‚
0.1 J/(cm2 · sec · K), X0 = 4.5 · 10-6 mole/cm3 [in-
ôÅ‚
Bq(-"H)k(T )X0
ôÅ‚
òÅ‚
Á(T ) = side the limiting cycle, there is one unstable steady
cp(q + V k(T ))(k(T ) + 2q/V )
state (thick curve); the arrows show the direction of
ôÅ‚
ôÅ‚
motion along the phase trajectories].
ôÅ‚
hS
ôÅ‚
ôÅ‚
- .
ół
cp(2q + V k(T ))
2
Here A = hS/cpÁ, B = E/RT , E is the activation passing through the boundary, the steady state changes
energy, and R is the universal gas constant. its type of stability. Aligning the curves L" and LÃ on
Figure 2 shows some bifurcation curves obtained by the same parameter plane, one obtains a so-called para-
formulas (8) and (9). Several multiplicity curves L" in metric portrait of the system. The relative position of
the (X0, V ) plane obtained by varying a third param- the curves L" and LÃ (see Fig. 3) determines the divi-
eter (T0) are shown in Fig. 2a. There are three steady sion of this region into subregions with different number
states inside the region bounded by the curve L", and of steady states (one or three) and different type of their
there is one steady state outside this region. One can stability (stable or unstable steady state). In Fig. 3, one
see from Fig. 2 that the region of multiple steady states can see the region with a unique and unstable steady
has the largest dimensions for relatively low inlet tem- state (inside the loop of the neutrality curve LÃ), which
peratures, all other factors being equal. The dimensions is responsible for self-oscillatory solutions of system (1).
of the region bounded by the neutrality curve LÃ (see Outside the regions bounded by the curves L" and LÃ,
Fig. 2b) also depend strongly on the parameter T0. In the steady state is unique and stable.
638 Bykov and Tsybenova
Each subregion of the parametric portrait of system
(1) (see, for example, Fig. 3) has its type of phase por-
trait which reflects the specific features of the dynamics
of the system and the relative position of its phase tra-
jectories for varied initial data. Figure 4 shows a phase
portrait corresponding to the case of unique and unsta-
ble steady state, which ensures the existence of a lim-
iting cycle for system (1). If the initial data lie outside
this cycle, all phase trajectories are  wound on it as
t ".
If it is necessary to construct bifurcation curves of
the type (8) and (9) in other parameter planes, one can
use the geometrical procedure proposed in [1, 2], which
is based on explicit parametric relations of the type
(5) and (6), or the numerical parameter-continuation
Fig. 5. Parametric portrait in the (X0, V ) plane for
method, which we modified for the case of one equation
model (1), constructed from the parametric portrait
with a parameter G(x, Ä…) = 0.
for model (10) for h = 0.055 J/(cm2 · sec · K) (the bi-
furcation curves L" and LÃ are similar to the curves
in Fig. 3).
RELATION BETWEEN DIMENSIONLESS
AND DIMENSIONAL MODELS
we obtain the bifurcation values of the dimensional pa-
With the parametric portrait of a model with di-
rameters V and X0. For these values, we construct the
mensionless parameters, one can readily construct the
parametric portrait in the (V , X0) plane.
parametric portrait in the plane of dimensional param-
Figure 5 shows the parametric portrait for
eters. For illustration, we consider the Aris Amundson
model (1) constructed in the (V , X0) plane by ana-
model widely used in the theory of chemical reac-
lyzing the corresponding dimensionless model (10) in
tors [5, 6]:
the (Da, ²) plane. In contrast to Fig. 3, the paramet-
ric portrait in Fig. 5 is constructed for slightly different
dx
= f(y)(1 - x) - x,
values of the parameters. In this case, the region of
dt
(10)
self-oscillatory regimes (inside the loop of the curve LÃ)
dy is large; the curves L" and LÃ are constructed with
= ²f(y)(1 - x) - s(y - 1).
the use of the parametric portrait of system (10) given
dt
in [7, 8]. The parametric portrait can be constructed
Here x and y are the dimensionless concentration
in other planes of dimensional parameters similarly; in
and temperature, respectively, and f(y) = Da exp(Å‚(1
some cases, however, technical difficulties arise due to
- 1/y)). Introducing the dimensionless parameters
the fact that the same dimensional parameter can en-
E hS
ter several dimensionless parameters. In this case, it is
Å‚ = , s = 1 + ,
RT0 qÁcp
necessary to perform a direct parametric analysis of the
dimensional model (1).
V (-"H)X0
Da = k(T0), ² = ,
q ÁcpT0
one can obtain, for example, DETERMINATION OF
IGNITION BOUNDARIES
q ÁcpT0
V = Da, X0 = ². (11)
k(T0) (-"H)
To ensure the technological safety of chemical pro-
To construct the parametric portrait in the (V , X0) cesses, it is necessary to determine, among other factors,
plane with the use of the parametric portrait in the di- the limits of existence of low- and high temperature
mensional (Da, ²) plane and formulas (11), it is neces- regimes. An analysis of time dependences and corre-
sary to specify the dimensionless parameters Å‚ and s sponding phase portraits shows that a steady state can
for which the last portrait is obtained. Substituting be attained with a considerable dynamic  throw. Even
the calculated values of Da and ² into the relation be- though a low-temperature steady state exists, there is
tween dimensionless and dimensional parameters (11), a wide range of the initial data for which the transi-
Parametric Analysis of the Models of a Stirred Tank Reactor and a Tube Reactor 639
tion of the system to the steady states is accompanied
by intense heating during the reaction. Therefore, we
consider dimensional and dimensionless models to find
a region of technologically safe regimes characterized by
monotonic transition to steady states in one of the pa-
rameter planes.
We consider the Zel dovich Semenov model for a
first-order reaction [1, 2]. The corresponding paramet-
ric portrait is shown in Fig. 6 in [1]. An analysis of
the phase portrait of the dynamic system with a low-
temperature steady state shows that there exists a re-
gion of the initial data for which the transition to the
steady state is monotonic. Outside this region, the tran-
sition of the system from the initial data to the steady
state is accompanied by considerable heating. It is quite
possible that abrupt heating of the system can lead
to an explosion or fire. Therefore, the line that sep-
arates regimes with low- and high-temperature steady
state may be called the boundary of technological or fire
safety.
For illustration, we also calculated the dimensional
model of a continuous stirred reactor (1) in which an
exothermic reaction occurs. Figure 5 shows the corre-
sponding parametric portrait in the (V , X0) plane. The
technological-safety boundary in the chosen parameter
plane is similar to that shown in Fig. 6 in [1]. It cor-
responds to relatively low values of the inlet concentra-
tion X0 and volume V , i.e., it lies below and to the left
Fig. 6. Temperature distribution (a) and concen-
of the curves L" and LÃ. Since the bifurcation curves
tration distribution (b) along the reactor [model
can be constructed by explicit formulas of the type (8)
(12)] for u = 9.171975 cm/sec, T0 = 583 K, and
and (9), one can study the effect of various characteris-
X0 = 4.5 · 10-6 mole/cm3.
tics of the reactor on the position of this curve.
Thus, an analysis of the parametric and phase por-
traits of the corresponding mathematical model allows
In the stationary case, the mathematical model of a
one to determine regimes of technological safety of the
continuous tube reactor in which an exothermic reaction
processes.
occurs has the form
dX
u = -k(T )X,
CONTINUOUS TUBE REACTOR
dl
(12)
In addition to the above-considered model of a dT 4h
cpÁu = (-"H)k(T )X + (Tw - T ),
continuous stirred reactor, the model of a continuous
dl d
tube reactor describing processes in tube reactors is
where l is the current length of the reactor, u is the
widely used in mathematical modeling of combustion
velocity of injection of the reactive mixture, and d is
processes. It is, therefore, of interest to compare the
the diameter of the reactor.
dynamic characteristics of a continuous stirred reactor
In model (12), the current length of the reactor
and continuous tube reactor in which the same reaction
varies within the limits 0 l L, where L is the total
occurs, provided the thermal and geometrical character-
length of the tube reactor. The inlet conditions have
istics of the reactors are similar. As an example, we con-
the form
sider the mathematical models for a first-order reaction
X(0) = X0, T (0) = T0.
and find the relation between the dynamic characteris-
tics of a continuous stirred reactor and the concentra- The volume flow rate is related to the velocity of the
tion and temperature profiles in a stationary continuous mixture by the formula u = 4q/Ä„d2. The values of
tube reactor. thermal and physical parameters and the ranges of their
640 Bykov and Tsybenova
variation are taken to be the same as for the model of edge of the dynamic processes occurring in a continu-
a continuous stirred tank reactor (see Table 1). ous stirred tank reactor allows one to estimate the tem-
The diameter and length of the tube reactor are perature profiles and critical conditions in a tube reac-
chosen from the condition that its volume is equal to the tor. The position of the hot spot and the magnitude
volume of the continuous stirred tank reactor. Figure 6 of heating in the tube reactor correlate well with the
compares the self-oscillatory regime in the continuous amplitude and period of oscillation of the temperature
stirred tank reactor with the corresponding tempera- in the corresponding spherical continuous stirred reac-
ture and concentration profiles in the tube reactor. The tor. Consequently, the laboratory data on the dynamics
self-excited oscillations in the continuous stirred reactor of combustion processes, which are usually obtained in
are characterized by periodic increase and decrease in the Frank-Kamenetskii well-stirred reactor, can be used
temperature, and the corresponding temperature pro- for approximate estimation of the characteristics of the
file in the tube reactor has a pronounced hot spot. The same process in the tube reactor.
calculations show that the position of the hot spot can
be determined approximately by the formula l" = ut",
where t" is the period of self-excited oscillations in the
REFERENCES
stirred reactor. Thus, the knowledge of the dynamic
characteristics of the stirred reactor allows one to esti- 1. V. I. Bykov and S. B. Tsybenova,  Parametric analysis of
mate the position of the hot spot in the tube reactor. the simplest model of the theory of thermal explosion 
Certainly, this estimate is rather approximate, since it the Zel dovich Semenov model, Fiz. Goreniya Vzryva,
depends on the type of the reactor and kinetic specific 37, No. 5, 36 48 (2001).
2. S. B. Tsybenova,  Parametric analysis of the basic mod-
features of the process. However, our calculations show
that the quantity l" may be used as a first approxima- els of the theory of chemical reactors and combustion the-
ory, Candidate Dissertation in Tech. Sci., Krasnoyarsk
tion to determine the position of the hot spot.
(1999).
It is noteworthy that the temperature of the hot
3. V. I. Bykov, T. P. Pushkareva, and Ya. Yu. Stepan-
spot is highly sensitive to variation in the parame-
skii,  Simulating self-excited oscillations in the cold-flame
ters [5]. There exists a peculiar bifurcation value of
combustion of a mixture of n-heptane and isooctane
parameters for which a smooth temperature profile be-
in a well-stirred reactor, Fiz. Goreniya Vzryva, 25,
comes a profile with a sharp rise in temperature, which
No. 2, 21 27 (1987).
can be interpreted as a thermal explosion in the tube.
4. V. I. Bykov and T. P. Pushkareva,  Simulation combus-
These profiles are shown in Fig. 6, where the coefficient
tion of a mixture of n-heptane and isooctane in a cylin-
of heat transfer through the reactor wall is the varied
drical reactor, Fiz. Goreniya Vzryva, 26, No. 2, 34 37
parameter. A similar effect is produced by varying the
(1989).
inlet concentration of the reacting substance.
5. R. Aris, Introduction to the Analysis of Chemical Reac-
tors, Prentice-Hall, Englewood Cliffs (1966).
CONCLUSIONS
6. M. Holodniok et al., Metody Analyzy Nelinearnich Dy-
namickych Modelu, Praha (1989).
Mathematical models of specific physicochemical
7. V. I. Bykov and S. B. Tsybenova,  Parametric anal-
processes contain dimensional variables and parame-
ysis of some basic models of the combustion theory,
ters that characterize a given geometry of the reactor
in: Chemical Physics of Combustion (Collected scientific
and particular physical and kinetic properties of the
papers dedicated to 70th Anniversary of Academician
reaction. In this case, a parametric analysis of the
G. I. Ksandopulo) [in Russian], Almaty (1999), pp. 133
dimensional model can be performed directly for var-
135.
ious combinations of dimensional parameters or indi- 8. S. B. Tsybenova,  The Aris Amundson mathematical
rectly by using the results of the parametric analysis of
model and its parametric analysis, Vestn. KGTU (col-
dimensionless models. For example, at least the fol- lected scientific papers), Krasnoyarsk (1970), pp. 137
lowing three dimensionless models correspond to the
139.
basic model of the combustion theory (1), which de- 9. B. V. Vol ter and I. E. Sal nikov, Stability of Modes of
scribes one exothermic reaction in a continuous stirred Operation of Chemical Reactors [in Russian], Khimiya,
tank reactor: Zel dovich Semenov model [1, 2], Aris Moscow (1981).
Amundson model [5, 6], and Vol ter Sal nikov model [9].
The results of the bifurcation analysis of these mod-
els can be used to estimate the critical conditions in
the initial dimensional model (1). Moreover, the knowl-