The Chemical Engineering Journal, @ Elsevier Sequoia S.A., Lausanne.
9 (1975) 215-222 Printed in the Netherlands
Use of the Heaviside Function as a Substitute for Arrhenius’ Law in Mathematical Models of Distributed Parameter Reactors Part II: First Order Reaction GUILLERMO Departamento (Received
E. HOUGH and ITALO H. FARINA de Ingenieria Quimica, Facultad de Ingenieria.
3 June 1974; in revised form
Universidad National de VIAPlota, 1 esq. 47-La
Plats (Argentina)
10 March 1975)
(0) Mass balance at the reactor wall:
Abstract The analysis of a previous paper is extended to the more interesting case of a first order reaction. The results are similar to those of the previous work. Partially ignited solutions are also considered. The assumption of a non-conductive wall needs to be revised in order to fit the model to actual reactors.
dCw Pz=-rtk,(Cf-G)=O
(1)
This balance indicates that the number of moles of A1 that disappear through chemical reaction is equal to the number of moles of Al that are transferred from the fluid to the reactor wall.
(6) Mass balance in the fluid phase: 1. INTRODUCTION
The present work is based on a mathematical model for fixed-bed reactors originally proposed by Aris’l * and later modified by Farinas. We shall refer to this last work as Part I. In this work our objectives, and the reasons for pursuing them, are the same as in Part I but here the comparison between the use of the Heaviside function and Arrhenius’ exponential is extended to first order chemical reactions.
A F
(2)
r = kog(W’)G Combining
(3)
eqns. (l), (2) and (3) we obtain
‘3 =- V-4 ‘3 - ko f(W’)Cf
(4)
where
We shall make the following assumptions. (1) The reaction takes place at the wall of the tube; this implies no mass accumulation on the solid phase. (2) There is a first order irreversible chemical reaction (Ar + Az). (3) The fluid is in piston flow. (4) There is no axial heat conduction and no axial mass dispersion. With these assumptions are as follows.
- k,P(Cf - C,)
Bearing in mind the supposed reaction order, the expression for r will be
A 2. MODEL
= - VA $
the mass and heat balances
d W’>
f(w’)=(k,,/k,)g(W’)
+1
(5)
(c) Heat balance at the reactor wall: cwAw $’
= Ph(T’ - W’) t (-AH)Pr
(6)
(d) Heat balance in the fluid phase: cfA $
t VAcf $
= Ph(W’ - T’)
(7)
216
G. E. HOUGH, I. H. FARINA
These equations can be made dimensionless defining the following variables, as in Part I :
The boundary
by
hPz’
Consideration side
aT
=mw (1
-x)
(8)
aT (9)
ar+az=w-7 $=T-
0
we easily find that (14)
of eqn. (13) shows that the left-hand
QI = h's - Ts
where Q = Pkc(-AH)Cr. T; is the critical or reference temperature and its meaning has been analysed in part I. Introducing these variables into the balance equations we obtain
ax ax
the above equations
z =
Xs =p(Ts -~ Ti)
hPt’ *=-Acf
Ci - Cf X=----ci
jy+z
are
T,(O) = Ti and Xs(O) = 0 at Combining
z =-VAcf
conditions
W+f(W)(l
(15)
is a measure of heat dissipation, hand side QI, =f(K)(l
whereas the right-
~ XS)
(16)
is a measure of heat generation at the reactor wall. Before proceeding with our study of the steady state we shall list the possible forms of the f(W’) function. We can consider three cases. (a) f(W’) is defined by eqn. (5). (b) When the mass transfer step is very fast, i.e. k, > ko, eqn. (5) reduces to f(W’) =gP_(W’)
-,I’)
(17)
(10) fun’,‘t’,o~IW’) is app roximated
by the Heaviside step
The two parameters are &Z/!
w-
&(W’) -g,(W;nax)H(W’)
cwAw cfA
if W’ > Tf
fi is proportional
to the inverse of the adiabatic temperature rise and w relates the heat capacities of the solid and fluid phases. The initial and boundary conditions for eqns. (8)-( 10) are T(0, t) = Ti(t)
W(z, 0) = We(z)
X(0, t) = 0
T(z, 0) = To(z) X(z, 0) =o
3. STEADY
STATE
In the steady state, eqns. (8)-(10) pectively,
(18)
where w,,, is the value of W(0, a). The meaning of Wmax will be clear when we consider totally ignited solutions. We can see that the correct approximation involves multiplying the Heaviside function by ge( W&,), which is the maximum heat generation obtainable. We shall now analyse the steady state for these three cases. ge( W’) given by Arrhenius’law Using the dimensionless form and defining the parameters
become, res-
dX, x = Pf(WSX1 - XS)
(11)
dT, = W, ~ TS dz
(12)
WS- TS =f(Ws)(l
if W’GTA
EhP (y=_
RQ
g,(T;)
cl=(ko/k)g,(T:)
t 1
eqn. (5) becomes
~
Xs)
,,=[z+exp(( (13)
z+$)-‘)]-’
(19)
THE HEAVISIDE
FUNCTION
AS A SUBSTITUTE
FOR
ARRHENIUS’
where E,
=
1-
elko/k,
(20)
El
Returning to the problem of finding a solution of the steady state given by eqns. (1 1)-(13), we see that an analytical solution is not possible because of the presence of the non-linearity given byf(W). Notwithstanding, we can reach some conclusions from the study of the heat generation and heat dissipation functions. These are coupled by eqn. (14). We shall first consider some of the restrictions that apply. The physical condition that W’ and T’ can never take values below 0 “K can be expressed by W and T > -o/in
E’
(21)
Another restriction arises from eqn. (14) and the fact that 0
(22)
In Part I it was seen that for (Y> 0.544 there cannot be three intersections of QI with QII for f(W) =g,(!Q. A consideration off(W) leads again to the limit for (Y;for if K = ko/kg + 0, f(W) +g,(W) and an increase in the value of K “flattens out” the f(W) curve. In order to study the more interesting case in which multiple steady states occur, we chose the following values of the parameters: CY= 0.1, ~1 = 0.3, K = 0.1 and /3 = 10. The values of ~1, K and fi were chosen to give a good illustration. The QI and QII curves for the above values of the parameters are
Fig. 1. Heat generation a=O.l,el =0.3,K=O.l
and heat dissipation andP=lO.
functions
for
217
LAW. II
represented in Fig. 1. It can be seen in Fig. 1 that for a given value of X point A, where the QI and QII curves cross, represents the steady state of the position in the reactor where the conditions of this point obtain. Since A is an arbitrary point, we can define a zone in the plane where we can find any one of the infinite possible steady states. That zone has been shaded in Fig. 1. Very rapid mass transfer This supposition leads to f(W) = g,(W), as was seen above. Defining the parameter e2 =g,(T:) we obtain 1 -‘I W-01 In.52 ii 1 I
g&V)=exp
I
(23)
Here ~2 has the same meaning as the E in Part I. Hence QII =g,(W)(l
- xs)
(24)
Thus the plot of Qrr with Qr will have the same form as Fig. 1, and the conclusions to be reached will be identical. Heaviside step function If we write eqn. (18) as a function variables we obtain
ge(Wmax) ifW>O &%(W)=
Fig. 2. Steady p= 10.
0
state solution
ifW
usinggh(W)
of the dimensionless
(25)
for fi = -10.08 and
218
G. E. HOUGH.
If we now introduce gh(W), for W > 0 (ignited reactor), in the steady state equations we can obtain an analytical solution of the steady state because the non-linear factor has been eliminated. Solving the resultant system gives
Xs = 1 - ew(-bd~,,,)~zl
(26)
T,l=Ti.+ (llP>[l - exp{-g~(~~~~)~z~l
(27)
WS = ri + l/S + exp C-ge(Wmax)PZl
ke(Wma~>
-
l/PI
(28) We can now plot these values of IV,, TS and X, (see Fig. 2). Later they will be compared with the numerical solutions for the transient state.
4. TRANSIENT
STATE
I. H. FARINA
w 2
1
t-0
o
Tt
O’tYt-
OF‘-
Unless otherwise indicated, the transient state will be analysed for constant initial and boundary conditions, which will be of the form T(Z,
0) = T(0, 2‘)= Ti
W(z,O)=
wo
(29)
X(0, t) = X(z, 0) = 0
Fig. 3. Temperature and Wo=O.
and concentration
profiles
for Ti = - 0.01
profiles
for Ti = 0
The values adopted for the parameters will not be changed: ~~=O.l,fl=lO,e2 =0.3ando= 1. 4. I. Analysis of the transient state considering mass transfer to be very rapid We can obtain different solutions depending on the values of Wo and ri. 4.1.1. Totally ignited solutions We consider a solution to be totally ignited when X> 0.1 for all z > 0.05 once the steady state has been reached. It should be noted that this definition and the ones to be given for quenched and partially ignited solutions are not completely arbitrary since they are quantitative expressions of the meaning contained in the respective name. We shall now analyse the profiles obtained. The profiles for rt = -0.01 and Wc = 0 are plotted in Fig. 3. With regard to these we observe the following. The maximum value of W is at the reactor entrance since it is here that the greatest heat generation through chemical reaction occurs. The value of W decreases because less and less heat is generated as the conversion progresses; it becomes practically zero
t
‘IWz 1
Fig. 4. Temperature and IV0 = -0.08.
2
3
and concentration
4
THE HEAVISIDE
FUNCTION
AS A SUBSTITUTE
FOR
when the value of X tends to one. Thus, since no heat is generated, W tends to equal Tat the reactor outlet. The overshoot observed in the T profile for r = 1 results from the fluid temperature being increased by heat transfer from the wall. Since the value of W is greater at the reactor entrance, the dynamics of heat transfer provoke the overshoot, although a continuous rise is observed in the steady state. The profiles for Ti = 0 and We = -0.08 are plotted in Fig. 4. The interesting point to be observed in this case is that with a very low initial wall temperature the reactor is ignited due to the “hot” inflowing current. 4.1.2. Quenched solutions We consider a solution to be quenched when X is less than 0.1 for all values of z after the steady state has been reached. This type of solution is a result of taking very low values of Ti and We, e.g. Ti = W0 = = -0.08. The heat generation is then very low and does not allow a sufficient rise in temperature to ignite the reactor. 4.1.3. Partially ignited solutions We consider a solution to be partially ignited when x<
0.1
219
LAW. II
w 1
t
-1t
0.1 t
-0.1
t
-co 25
Fig. 5. Temperature and W0 = -0.07.
and concentration
profiles
for ri = -0.07
forO
and x>o.9
ARRHENIUS’
for ze
where 0.05
As the time is increased the conditions are attained for lower values of z, and ignition occurs at z = 0.8 for the steady state. The profiles stop at this point because the fluid needs a minimum time (that necessary to cover the distance z = 0.8) to reach the ignition temperature. For the totally ignited solution this time is zero because the temperature at the reactor entrance is sufficiently high to allow ignition. We observe a discontinuity in the profiles. This is because axial heat conduction in the solid phase has not been considered.
1
a9
x-
E-
0
0.5
1
Fig. 6. Heat generation
w
functionge(W).
220
G. E. HOUGH,
4.2. Profile movement An important feature of the partially ignited solutions is the movement of the profiles after they acquire the steady state form. This movement must not be confused with the “creeping profile” concept described in the literature, which refers to the movement of the profiles, at a regular speed, when the steady state is perturbed. This feature is not expected in our model because of the lack of diffusive terms, as is explained in the following. Later we shall see that the shapes of the steady state profiles using either gh(W) or ge(W) are very similar. We have also seen that the use ofgh(W) allows an analytical solution of the steady state, given by eqns. (26)-(28), to be obtained. These are the equations we shall now use to estab!ish the initial profiles in the reactor, for ignition at an intermediate point. Figure 7 shows the X profile (T and W behave in
‘:_&_@_
I. H. FARINA
____
I
:
0°K -007
Fig. 8. Quenched
-0 05
-Tt
and ignited
zones
ignited solutions. The rest of the plane corresponds to totally ignited solutions. Zone A was obtained by observing the steady state reached through solving the transient state for various values of r, and WO. Zone B was determined by solving eqn. (10) at the reactor entrance, as will be seen in Section 4.4. It can be seen from Fig. 8 that both the A and B zones are relatively small. This can be predicted from Fig. 6, where we see that the zone in which three intersections (one of which is quenched) between Qt and Qrr are possible is very small.
,rZ
1
Fig. 7. Concentration in Ti.
3
4.4. Wall temperature at the reactor entrance If we consider eqn. (10) at the reactor entrance, taking X = 0 and Ti = constant, we have
4 0
profile
movement
after a perturbation
an identical way) for t = 0, with ignition occurring at zc = 2.5. We can see how the profiles evolve with time for an inflowing current at T(0, t) = -0.065. Since no diffusional effects exist, the reactor zone between z = 0 and z = zc acts independently of whether the rest of the reactor is ignited or not. This explains why the profiles develop in the way described for partially ignited solutions. Hence the profile for t = 0 moves to another zone for the steady state, but not through a creeping effect. For a cold inflowing current (T(0, t) = -0.08) no variation from the t = 0 profile is observed. This is because the fluid does not cool the reactor wall sufficiently, in the zone where heat is generated, to displace the profile towards the reactor outlet. 4.3. Ignition and quenched zones in the WO, Ti plane In a plot of Ti versus WO (Fig. 8) different zones can be identified according to the type of solution obtained from the values of T, and WO within these zones. This plot corresponds to Fig. 6 of Part I. Zone A indicates quenched solutions and zone B partially
dWi o_dt=
WitTiff
Wi(0) = WOz
(30)
The solution of this equation was obtained for various values of Tt and WoZ and the results are plotted in Fig. 9. Here we observe that, for Ti = constant, the
b /.
Fig. 9. Wall temperature
0 08
at the reactor
entrance
THE HEAVISIDE
FUNCTION
AS A SUBSTITUTE
FOR
value of Wi follows a vertical line until steady state is reached. There are two curves for the steady state:
t I=- t F
ARRHENIUS’
w
221
LAW. II
=m
curve (a), corresponding curve (b), corresponding
J3
to an ignited state, and to a quenched state.
0 When the solution evolves to curve (a) we can be sure of obtaining a totally ignited solution for the reactor. When the result is curve (b) (which is the line Wi = T$ the solution for the reactor can be quenched or partially ignited. With this in mind it is easy to see how zone B of Fig. 8 was determined through the wall temperature at the reactor entrance. In Fig. 9, for a given Wo, the arrow indicates the minimum value of Ti for which the solution evolves towards curve (a). 4.5. Analysis of the transient state using the Heaviside finction We shall consider two cases involving different values of Ti and Wo. The profiles for Ti = -0.01 and Wo = 0.01 are plotted in Fig. 10. This case is comparable with the first case in Section 4.1.1. We can see that for WO = 0 use ofg,,(IV) will lead to a quenched solution. Comparing the profiles in Figs. 3 and 10 we see that the differences in both the transient and steady states are very small.
1 -1
1 0.1
1
2
3
2,4
=
t
0
=a
Jo
1
2
3
2
4
-5 13
,.I/
Fig. 11. Temperature G = 0.01, Wo = -0.08
and concentration profiles and ge(Wmax) = 0.9046.
for
The profiles for Ti = 0.01 and WO = -0.08 are plotted in Fig. 11. This case is comparable with the second case in Section 4.1 .l . From Fig. 4 we can see that for t = 1 there are values of X greater than zero. In Fig. 11, however, it is clear that values of X greater than zero only appear when the inflowing current heats the reactor wall above the value zero. Likewise the discontinuity which is observed in the X profile for t = 3 can be attributed essentially to the discontinuity present in the step function. We see that gt,(W) gives a good approximation for the steady state. Both these cases are totally ignited solutions. Quenched solutions are obtained by taking both WO and Tr to be less than or equal to zero. Partially ignited solutions cannot be obtained using gh(W) with constant boundary and initial conditions.
5. CONCLUSIONS
Fig. 10. Temperature and concentration Ti = - 0.01, Wo = 0.01 and g,(Wm&
profiles
= 0.9026.
for
With the introduction to Part I in mind, it is of importance to analyse to what degree the Heaviside step function is a good approximation to Arrhenius’ law. The greatest advantage, of course, is the computer time saved by using the former. For our case the
G. E. HOUGH,
222 computer time was reduced by approximately onehalf when ge(w) was replaced bygh(W). We shall now make a comparative analysis of the results obtained using&,(W) andg,(W). (1) Firstly, as already mentioned in Section 4.5, partially ignited solutions cannot be obtained using gh(W) with constant initial and boundary conditions. (2) For Wc > 0 both the transient and steady states obtained withgh(W) are good approximations compared with those obtained usingg,(W). (3) For We < 0 and z > 0 we have a good approximation in the steady state, although this conclusion is not valid for the transient state. (4) Perhaps the greatest difference in usinggh(W) is found in the Wc and Ti values that ignite the reactor. As we saw in Fig. 8, the quenched solution zone is very small if we use g,(W), but it is greatly increased whengh(W) is used, occupying the whole of the third quadrant in the Wc, F, plane. As a general conclusion from the above comparisons, it is clear that by using gh(W) a good approximation is obtained for the ignited steady state. For the transient state and other types of solutions the discrepancies indicate that the practical use of the Heaviside step function is limited. Two main problems remain to be solved for this model. In the first place an adequate adjustment of the parameter values would make possible a better representation of real reactors with this model. For example, we think it would be of interest to find the relation between Aw and the reaction area of a real fixed-bed reactor. Secondly, the inclusion of diffusive terms would lead to a marked change in the behaviour of this reactor model, mainly with regard to the discontinuities observed in the wall temperature profiles. Furthermore, the inclusion of axial diffusion would change not only the shape of the solution but also the number of steady state solutions.
I. H. FARINA
ACKNOWLEDGMENT
The support given to G.E.H. by the Consejo National de Investigaciones Cientificas y Tecnicas of Argentina, which made this research possible, is greatly appreciated.
NOMENCLATURE
A
C C AH h kg ko P T t W X z
cross section concentration volumetric heat capacity heat of reaction heat transfer coefficient mass transfer coefficient pre-exponential factor of Arrhenius’ law internal reactor perimeter fluid temperature time coordinate reactor wall temperature conversion length coordinate
Indices
f i S W
fluid entrance steady state wall dimensioned
physical value
REFERENCES 1 R. Ark and 0. Schruben, Chem. Eng. J., 2 (1971) 179 2 I. H. Farina and R. Ark, Chem. Eng. J., 4 (1972) 149. 3 I. H. Farina, Chem. Eng. J., 9 (1975) 207.