- Email: [email protected]
Theory of heterogeneous combustion instabilities of spherical particles
THEORY OF HETEROGENEOUS COMBUSTION INSTABILITIES OF SPHERICAL PARTICLES NORBERT PETERS Institut f~r Thermo- and Flulddynam~k, Techn~sche Unlversit~t Berlin, Germany The linear and nonlinear stability characteristics of the heterogeneous combustion of spherical particles are investigated on the basis of a simplified mathematical approach using integral relations. A condition for instability is derived which relate~ the parameters of the problem in an algebraic inequality and reflects the inttuence of internal diffusion and reaction. In a ease where three steady states exist only the lower one was found to be stable to infinitely small and to finite disturbances. Calculations of the transient behavior of the combustion of carbon particles are able to explain the nature of experimentally observed oscillatory instabilities. They appear to be caused by the unsteady heat exchange between the surface and the interior of the particle which produces a time lag. At large values of the thermal conductivity inside the particle the oscillations are damped and stability is obtained. that in the combustion of non-porous carbon particles mlstable steady states may exist, even The ability of heterogeneous combustion proc- if heat conductivity inside the particle is assumed esses to generate combustion instabihties such as to be vanishing. The question arises whether the internal difquenching, ignition and oscillations is known to be due to the interaction of heat release by chemi- fusion is an essential instability factor. For incal reactions with heat conduction and radia- stance, the recently observed dynamic combustion. If the stationary heat balance at the surface tion instabilities 15 of dense carbon spheres could is plotted against the surface temperature,~ not fully be explained by the influence of internal it can be shown that one or more steady state diffusion. In another case, where oscillatory h~solutions of the problem exist, some of which stabilities of the heterogeneous hydrogen-oxygen reaction on a pl~thlized spherical pellet were obmay be unstable. In recent years much effort has been concen- served,16 the role of internal diffusion could be trated on investigating the stability of single excluded as the cause of these instabilities. 9 catalytic particlesY -~3 Complex mathematical W]dle most of the above theoretical works promodels have been elaborated in these studies to vide considerable understanding of the instabilipredict stability aud uonstationary heat and mass ties of heterogeneous reactions on spherical partitransfer inside the particle, while transport from cles, the highly sophisticated nature of these the outside to the surface was imposed by assum- solutions makes their application by most ening constant Nusselt and Sherwood numbers in gineers rather improbable. The theory presented the surrounding medium. here is intended to provide a simple analytical In heterogeneous combustion processes, how- tool to explain heterogeneous combustion inever, nonstationary heat conduction, radiation stabilities. The mathematical description is simand diffusion in the surrounding gas phase may plified by using integral relations to reduce the have an important influence on the stabihty of system of partial differential equations to ordithe heterogeneous reaction. Thus it seems that nary ones. Integral methods, which have been the interaction of the reaction with the non- extensively used in fluid dynamics,~7 are very stationary transport processes in bottt phases helpful in problems where only information at the should be considered. A linear stability analysis boundaries of the region of integration is reby Wiuegardner trmt Sehmitz~ indicates, in fact, qub'cd. When used in stability theory, however, Introduction
363
364
HETEROGENEOUS COMBUSTION
~ i AH'R(T,c)
the most genera[ case, it is assumed to diffuse through the pores into the interior of the sphere where it reacts internally with the solid combustible material. All of these processes are essentiaRy time-dependent and will be investigated on the basis of the balance equations in spherically symmetrical geometry. For analytic convenience, the density, thermal diffnsivity, conductivity, diffusion coefficient and heat of reaction are assumed to be constant. These restrictions could be relaxed if necessary. Governing Equations
= T..T. ~ Ts w a. X ~W
The non-dimenslonalized balance equations describing the system are: inside the particle ( 0 ~ O0
T~
/0~0 200\
1)
1
Oc /02c 2 0 c \ ~= ubl~+; g)
I1~ Z LU
(j
(1)
AH.r(O,c)
I I
es Cs + ~s
P(0, c)
;
in the surrounding atmosphere
ro
RADIUS [cm] I~G. 1. Schematic of the process of non-stationary heat and mass transfer in the combustion of a spherical particle. the application of arbitrarily chosen trim functions restricts the number of possible pertubations causing instability. Therefore the aim of the present theory is to detect causes of combustion instabilities, rather than to assure the stability of a state of a system.
(l~v< ~ )
00 --~020 200 0r Or?~ ~ 07
0
0
r
0c
(3)
_ /02c. 2 0 c \
These equations are subject to initial conditions 0(7, 0) = e0(,)
(5)
c(~, 0)= ~(n)
(6)
and to boundary conditions at the center of the sphere, the surface and at infinity: ~?=0:
The Physical Model A single spherical particle of combustible material is considered in an infinite stagnant oxidizing atmosphere (see Fig. 1). A heterogeneous reaction X I + solid--*X, that converts the oxidizing gaseous component Xt into the gaseous product X~ takes place at the surface and in the interior of the particle. The heat of reaction generated on the surface is transferred by heat conduction and radiation to the surroundings and into the interior of the spheres. The oxidizing component is supplied by diffusion from the surroundings. In
oola~=o,
(7)
octan=o
'r Q = -- k,(00/0,7) = -- (00/0'~)+ QR(0)-'2d-/" R(0, c) (S) ./= --pLcb(0c/0,7) = -- Le(0c/e,1)+ R(0, c) (9) y----~ o0 :
O= 1,
c= 1
(10)
THEORY OF HETEROGENEOUS COMBUSTION Equations (8) and (9) represent heterogeneous balance equations that relate the ~wo phases. The radiation at the surface QR(0), the reaction rate inside the particle 17(8, c) and the heterogeneous reaction rate R(O, c) at the ~rface are defined as
F = kb exp(--3b/~)'C ~
(12)
R = k exp(--fl/O~,), c. ~
(13)
where reverse reactions have been neglected. Consider now a finite disturbance 0(~, r)=O(V, r)--O~(~)
(14)
5(~, r ) = c(~, 7)-c,(~)
(1~)
of the steady state solution 0,(~/), c,(v) of the system. When they arc introduced into Eqs. (1)(10), due to their linearity Eqs. (3)-(7) and (10) remain unchanged, while Eqs. (t), (2), (g), and (9) take the fona~
a~ de _
/a~. 20~. ~H.(r-r~) /o~5 2 oe )
were chosen taking into account the fact t h a t the finite disturbance must decay at large distgnces fronl the sphere:
[~O')+dt(r)(n-- 1)
O(n, r
(11)
Q~=r
(r-r~)
+&(r)(n-1)~] 9 expl--'Kn-1)]
(18)
-re~(z)(V--1)~] 9 e x p [ - - 7 ( ~ - - l ) ]
(19)
5(~, r
The exponential decay parameter 3"= 10.O was chosen such t h a t good agreement with the linear stability conditions of Winegardner and Sclunitz ~4 was obtained in the case of their particular example. When the coefficicnts in Eqs. (16) (19) are eliminated (~e Appendix), integration of the finite disterbance Eqs. (3a) and (4a) in the form
/,
o~d,=j,f-/o~ 2o~d
(2o
/ O'~C 2~5\
(2,
-
~
(la)
(2a)
Q = - ~,( o~/ o, ) = - (o~/o~) +(QR--Q~,)--AH(R--R,)
365
(8a)
yields the integral relations that determine the solutions. The weighting function ~ was introduced here to obtain analytical integrals. When this is done inside the particle the values close to = 1 are emphasized whereas the values at, ~/---~O are neglected. For this reason the differential equations at n = 0 together with the condition of symmetry in the limit */--~0
J = -- pLeo( O~/On) = -- Le( O~/&l)+ (R - - R~)
OO/O'r= 3a(O2~/ar
(22)
(9a)
Os/a~= 3LeKaW a~~)
(23)
Since integra2 relations are to be used to study the behavior of the above system, appropriate trial functions are required to approximate the solution in the domain of integration. In the present case it was found t h a t a Pohlhausen-type polynomial trial function could approximate the finite disturbance profiles inside the particle:
were used tot the integral relations inside the particle. When this is done there results a system of two second-order equations describing the time dependency of the surface temperature and concentration, and two first-order equations for heat transfer and diffusion flux at the surface
0(~, r ) = ~.(r)q-aKr) (1-- r/)+o~(v) (l--v) ~
BI(O2~da~)+ B2(Og~/O~)+ ~(a/o~)E(Q~ - Q..)
+aK'r)(1--n)~+a~(r)(1--~) '
(16)
e(n, r) = 5~(r)-4-bKr) ( 1 - ~/)+b~(r)(1--n) 2 -I-bdr)(1--n)~+b~(r)(1--n) ~ (17) In the surrounding gas phase the trial functions
- ~.
( R - R.)~-- B~(aQ/O~-) = Q - ~
-E(Qn-QRs)-Att.(R-R.)]
(24)
Bt(O~5~/Or:)+LeB~(c3~/Or)-l-Ba(O/Or)['R--R,] --Ba(OJ/OT)=Le(J'--Le'5.,--(R-- R.)) (25)
HETEROGENEOUS COMBUSTION
366
(7a/k,) (O0./Or)-'}-O~/ Oa--}-20a( O'O~/Ov)
(aH/t,c,)(o/o~.)(r- r,),=~
-
=
60(a~/l~)Q"F60(aAH/pc~,)G
(26)
(7/0) (OJ/O'r)-I-O~6~/0~"}- 20Loo(66~/Or) + (l/p) (0/0T) (F-- F~) ~-, =
60(Le~/o).]-6O(Leb/p)G
thus stability, while a positive real part of any one of the elgcnvalucs, signifies exponential growth of the disturbance and thus instability. The latter is the ease, as may readily be seen, if the coefficient p0 in Eq. (32) takes a negative value. From this condition, using the values of the matrix A, an algebraic criterion for instability may be derived after some elementary manipulations. The system is unstable if
(27)
AH. (OR,/OO+OGJOO)(,gR,/Oc+OGg'Oc) where
< {-- 1--OQ~cJOO+AH(ORjO0 G= K(r-
r~),=~+ ~ ( r - r , ) , ~
(28)
and B~= ( ~ + 3)/~ ~
B~= (3~+6)/v 2 B~= (3~-1-6~--6)/'y4 (29)
The second-order derivatives in Eqs. (26) and (27) may be eliminated by using Eqs. (24) and (25). When each of the second order equations is replaced by two first-order ones, tlmre results a system of six first order nonlinear differential equations t h a t are to be solved simultanconsly. The unknowns of this system may be written as a vector
9 = I0., ~, oOJo~, o~Io~, Q, J}
(3o)
which represents the disturbances of the steady state solution.
Linear Stability Analysis The stability of a steady state solution for infinitely small distances may be analyzed by linearizing the nonlinear terms in Eqs. (24)-(27) around a vanishing vector of disturbances. One obtains a sixth order matrix equation of the form
Ox/Or= Ax
(31)
The elemenk~ of the matrix A (see Appendix) are evaluated for the steady state. The eigenvalues of the system (31) may be obtained by solving its characteristic equation
hn-}-p~5"~-p4M-'~-psh3+p~?--F-pl;k.4-po=O (32) Negative real parts of the eigenvalues signify damping of an infinitely small disturbance and
~-OGdO0) }. { Le-}-ORJOc-I-OG,Oc}
(33)
It is interesting to note that neither the diffusion coefficient nor the thermal diffusivity or thermal conductivity inside the particle appears in Eq. (33). t51rthennore, by combining the reaction rates in the interior of the sphere with the heterogeneous reaction at the surface, an overJdl heterogeneous reaction rate may be introduced
Re*= R,+G,
(34)
If R, is replaced by R,* with G,= 0 in Eq. (33), the stability condition takes the same form as if only a heterogeneous reaction was considered. Thus the stability condition (33) remains unchanged if a combnstion model is adopted that neglects internal diffusion and assumes the interior reaction to takc place at the surface of the sphere by summarizing the reaction rates appropriately. Note t h a t Eq. (33) is a sufficient but not a necessary condition for instability to small disturbances. Using the Hurwitz stability criterion a complete ~ t of six inequalities may be derived. Numerical investigations show however that for a large range of parameters of interest, Eq. (33) is decisive for the change from stability to instability. As it is an algebraic relation which is easily applied to a given problem, it may be valuable for practical purposes.
Nonlinear Stability Analysis The above theory has been applied to the problem of combustion instabilities reported by Kurylko and Essenhigh. 15Iu their experiments on the combustion of a radiation-ignited spherical carbon particle these authors found oscillatory instabilities with temperature excursions between 25~ and 100~ to occur at furnace tempcratures of about 1000~ In view of the above results
THEORY OF HETEROGENEOUS COMBUSTION I/.00
i
r
t i iI
~L 1300
L/J
ik/ O_
367
1200
LU
i
]
I i '~
"
J
'i /
1100
','4-I
~A
~,o
1000 0,5
1.0
1.$
2.0
RADIUS [cm] Fro. 2. Steady state solutions. it was considered t h a t diffusion into the interior and internal reaction could be neglected in the analysis of this problem. Thus J and F are zero and Eq. (27) m a y be dropped in the system (24)(27). A system of five nonlinear first order differential equations was integrated numerically b y a Runge-Kutta procedure. Consistent initial values were obtained by determining the eigens vectors of the iinearized system (31). From the experimental data given~5 the following parameters could be established and used to calculate the steady state solutions (See Appendix) : r0 = 0.635 cm
T~ = 1000~
0~ = 1.22
Le = 1.5
f~ = 24.5
~,H= 10.0
The reaction order was assumed to be n ~ 1. A major difficulty consisted in choosing the preexponential factor ]c and the radiation coefficient (r. Figure 2 shows a plot of the steady state surface temperature over the radius r0 for three different values of k. The condition 0~= 1.22 at r0= 0.635 relates a value of g definitely to each value of ]z (see Table I). In Fig. 2 linearly unstable steady states are represented by a dashed
line. The largest eigenvalues for the three cases at r0= 0, 0,,= 1.22 are also given in Table I. I t is seen from Fig. 2 t h a t in case A, three steady states exist at ru=0.635, where only the lower one is linearly stable. In case B, only one steady state exists for all values of r0 which is linearly unstable at to= 0.635, while in case C, the singlc steady state at this radius is stable. A plot of steady state heat generation by reaction and heat transfer by conduction and radiation above the surface temperature is given in Fig. 3 for case A. The steady state solutions are those where the two liues cross e~ch other. By considering the slope of both lines, it is generally argued that the upper and the lower steady state is stable while the intermediate is unstable. Since linear stability analysis predicts instability TABLE I Case
k
~
;~(max)
A B C
3.4.10 ~ 1.0-10 TM 1.0.10 ~I
9.834 11.26 13.65
270.6 15.55-{-41.96 i -8.223
Evaluated at r0 ~ 0.635 cm, 9~ = 1.22
HETEROGENEOUS COMBUSTION
368 is
//
IO n., IM u, g,) z
r < uJ -t,..
~
O
R
§ QC
CASE A
/
o
-$
J
toflo
000
WALL
1200 TEMPERATURE
I1.00
[K]
FzG. 3, Steady state hea~ balsalce in case A.
1300-
ILl n,"
120@
ILl UJ I,--
1100-
100@
TIME
'r =
tag ro2
[_]
Fzo. 4. Transient behavior of the surface temperature s t ~ t i n g from the upper steady state sohztion in case A.
T H E O R Y OF H E T E R O G E N E O U S C O M B U B T I O N
369
y
IS
10
/ 1
AH- R
~/OR*
CASE B
J
Y
-w Boo
Oc
1200
tooo
1400
WALL TEMPERATURE [K] Fx~. 5. Steady state heat balance in case B.
50
^~AAA#
,,4 Z
AI
< < -SO 5O
2,0
2.5
3.0
~.0
L.S
$.0
2 <
o. X
w
~
-SO
3.S
TIME
I"-- ~
r2
5.5
[-]
FIO. 6. Transient oscillatory behavior of the surface temperature vadation in case B.
370
HETEROGENEOUS COMBUSTION 2S
----..
r162 .p uJ
-25
/
50 uJ
~X -50 0.0,~ ,
z 0
iJ 7 o
,,.
-o,0~
3.50
3.55
3.60
TIME v = t'ro~
3.65
3.70
[-]
FIG:. 7. Surface heat transfer, temperature and concentration variations of one oscillation in case B. for the upper steady state, a nonlinear calculation was carried out to determine its transient behavior. Figure 4 shows surface temperature against time as determined by the numerical solution of Eqs. (24)-(27). As the tcmpcrature disturbance was assumed to be only 0.03~ an induction period was required before the temperature dropped rapidly. The chemicaI reaction is quenched and aftcr an overshoot %o lower surface temperatures the stable steady state at 1042~ is reached at r = 3. This calculatioa suggests that the experimentally observed oscillations do not correspond to a situation where 3 steady states exist, since in this case the lower stabIe state would be reached. Furthermore it is shown that physical argu-
mentation may fail in certain cases while linearized mathematical analysis may predict tile correct physical behavior. If the steady state heat balance is considered for case B (Fig. 5), the slopes of the two curves differ considerably at the siugle steady state. Nevertheless the nonlinear calculations (Fig. 6) show exponentially growing oscillations after an induction period of r = 1.5 at ~he assumed initial temperatm'e pertnbation at 0.014~ At r = 4.0 the oscillations reach a constant limiting cycle with an upper amplitude of 52~ and a lower limit of 40~ To aoalyzc the oscillations, the cycle between r = 3.5 and r = 3.7 has been considcred in detail. Figure 7 shows the time-dependent behavior of heat trausfer, surface tem-
THEORY OF HETEROGENEOUS COMBUSTION
371
SO
" ~
T ~ , 3.62 "'
~ /
0
-so
.• o
T = 3.Se /,"-"--
t
2 RADIAL DISTANCE ;7 = r/r o
v
9 3,67
-SO
RADIAL DISTANCE -- r/r o
Fro. 8. SpaciM distribution of the temperature variation of one oscillation ia case B. perature and surface concentration. There is a phase shift between heat transfer to the surface and the su~ace temperature. Heat is continued to be furnished from the iuside of the particle at a nearly constant rate between r = 3.52 and v = 3.57 while tlie temperature drops continuously. Then, when the temperature begins to rise again, the heat transfer suddenly changes its sign. Heat is conducted at a high rate into the interior of the particle and thus stops the rise of the surface temperature. Prachcally no phase shift is found between the surface temperature and conceutratiou. These obscrvations lead to the conclusion that the experimentally observed oscillations are caused by thermal instabilities ratber than by internal diffusion. The temperature profiles inside and outside the particle for the same period of time are shown in Fig. 8. The amplitude maximum lies inside the particle at a certain distance from the surface. Heat exchange by conduction between this region
and the surface causes a time lag which is responsible for the instabilities. For the above calculation, values of ]c~=2.2 and a = 1.0 were asstoned. Since the value of thermal conductivity k~ is not needed to determine the steady state, it may be altered without changing the above considerations. When values of kt = 2.4 and higher were chosen, it was fmmd that the steady state was stable to finite disturbances. This indicates that the influence of porosity on the stability limits found by Kurylko and Essenhigh~b may be indirectly due r a change in thermal conductivity. This may explain the stability at low porosities. Case C was found to be stable to finite disturbances for all values of k,. From the comparison of the values of k and a iu the three cases considered, where only the intermediate one exhibits the possibility of oscillatory instabilities, it is concluded that there is only a narrow range of
HETEROGENEOUS
372
parameters in which the observed oscillations could occur. This is in agreement with the observations of Kurylko and Esseahigh ~ and may explain those limitations to the region of instability which were not investigated in the present analysis.
COMBUSTION
B r 7 9
6 Conclusions
1. Linearly unstable sLeady states may exist in the heterogeneous combustion on a spherical particle even if internal diffusion and combustion are neglected. 2. To detect the eau~s of these instabilities, nonstationary heat conduction and diffusion to the surroundings, i.e., non-stationary Nusselt and Sherwood munbers, must be considered. 3. If multiple steady states exist for a given problem, the assmnptiou t h a t the outer ones are always stable may not be valid. 4. Experimentally observed oscillatory ins~bilities can be explained b y a phase shift between the surface temperature and heat transfer from the interior. Nomenclature
a ab/% dimensionless thermal diffnsivity c cl/c1| mass fraction of component 1 cp c#/%g dimensionless specific heat E G
activation energy combined internal reaction rate Eq. (28) 511 (hl--l~)ct~/cm T~ dimensionless beat of reaction h enthalpy J j~ro/p#aocl~dimensionless diffusion flux k koaro(p~l~)"-l/aa dimensionless preexponential factor of the heterogeneous reaction kb "]~r dimensionless preexpcnential factor of the internal reaction kt ktb/kta dimensionless thermal conductivity Le D./a# gas phase Lewis number Le~ D~/% internal diffusion Lewis number n reaction order Q q~ro/k,~ T~ dimensionless heat transfer to the surface Qc 0 ~ - 1 dimensionless convective heat transfer Q~ dimensionless heat transfer by radiation Eq. (11) R a universal gas constant R dimensionless heterogeneous reaction rate Eq. (13) R* dimensionless overall heterogeneous reaction rate Eq. (34) r radial distance r0 radius of the sphere T Temperature t time
X
P o~t
r
E/RgT~ dimensionless activation energy dimensionless internal reaction rate Eq. (12) decay parameter emissivity for radiative heat loss r/ro dimensionless radial distance TIT** dimensionless temperature eigenvalue Pb/Pgdimensionless density ~teroT~/k,adimensionless radiation coefficient Boltzmann radiation constant agt/ro~ dimensionless time
Supexscript disturbance from the steady state
-
Subscript b q s w ~r
solid body surrounding ge~ steady state surface (wall) at infinity Appendix
The coefficients in Eqs. (16)-(19) are eliminated by using the boundary conditions at the origin and at the interface as well as by the requirement t h a t the trial functions mus~ fulfill the differential equations at 7 = 1. Finally, Eqs. (16) and (17) must be symmetrical around the origin, One obtains the following equations: al+2~+3az+4ad= 0
(All
l)1+2b2+ 3b3+ 4b4= 9
(A2)
~ = Q/k,
(A3)
b~= J/pLeb
(A4)
o~/o~=a(2ar 2a,)+(~/z~p)(r- r.)~.l (AS) o~/or= L~(2b2--2h)--(1/p)(r--r.)~.l (A6) a~-{-2a2-1-4 a s + 8a4= 0
h+~b~+4bs+Sbd= O
O=--(d,-7~)+(Q~--Qn~
(A7)
(AS)
R.) (AS)
] = -- Le(et-Te.)+ (R-- R,)
(A10)
O~/Or=26z+2d~(1--7)+~(7~--27)
(All)
0~/0~= Le[2~+2~(1--v)+~(v'--2x)] (A12)
T H E O R Y OF H E T E R O G E N E O U S COMBUSTION The elements of the matrix A, Eq. (31) have the following form: A(1,1)= A ( 1 , 2 ) = A ( 1 , 4 ) = A(1,5)= A(1,6)= 0.0
373
where
Bi = 60B{7,
B5 = 20B3/7,
B 6 = 20B1/7
A (2,1) = A (2,2) ~ A (2,3) ~ A (2,5) = A (2,6) = 0.0 A (3,6) = A (4,5) = A (5,6) = A (6,5)= 0.0
(A14) and
A(1,3) = A (2,4)= 1.0
D,=BI+k~Bj7a,
A (3,1)= ( l / D , ) { - - 1-- (0/oM)
D,=B1+pB~/7
(A15)
9[QR,--AH(R,+B4aG,)]}
The internal reaction rate F, in Eq. (A13) is evaluated at y = 1. The steady state solution for n = 1, F = 0 a n d J = 0 is given b y
A ( 3 , 2 ) = (1/D,) { AII (a/Oe)[ R,+ B4aG,]} A (3,3) = (1/Dr) { -- B2-- s3(a/ao)
9[QR,--AH(R,+F,/7)]--B~k,} A(3,4) =
O/D,)[Ba&H(a/ac)( R , + r J 7 ) ]
c,(~)=BJ~+l
A(3,5) = (1/Dt)(1--B4a)
(1_>,/>oo)
(A16)
where
a ( 4 , 1 ) = (1/D~) {-- (a/ao)
9(LeR~+B~LebG~) }
B~:
k- exp(-~/o~,) Le+k. exp(--fl/8~,)
(A17)
A ( 4 , 2 ) = (1/De)I-- Le~- (O/Oc)
O,(*/)= 0~,
9(LeR~+B4LebG~)] A (4,3) = (1/D~)[-- B3(O/00) (R,+ rdT)]
(0_>n_> 1)
(A18)
and
A (4,4) = (I/D,) [-- B2Le-- Ba(O/Oc) 0,=Bdn+l
9(R,+ r d 7 ) - B~pLeb] A ( 4 , 6 ) - (1/D~) ( Le-- B4Leb) A (5,1) = (1/D,)I (k,/7a)[l+ (a/ao)
( l > y > ~o)
where Bt= 0~-- 1
9 ( Q n , - &H R,)_-]+ B#AH(aG,/00) } A (5,2) = (1/Dr) I - - (kJ7a)AH(OR,/Oe) A(5,3) = (1/D~) {(k,/7a)[B~+B,(a/ao) x ( Q R , - A H R , ) ] + (B1/7)AH(OF,/OO)- Bdv,}
0~,-- l + a ( 0 ~ , ' - - 1)
AHLek exp(--fl/8~) ~
A(5,4) = ( l / D , ) [ - - (kd7a)BaAH(OR,/ae)
= o
(A21)
+ (B~/7)&H(aFo/Oc)] A (5,5) = (lIDs) (-- k~/Ta-- 3B~a)
REFERENCES
A(6,1) = (1/D~)[(pLe/7)(OR,/O0)
-- B4Leb(OG,/aO)] A(6,2) = (1/D~)I (p/7)[Le~+ Le(OR,/Oc)]
-- B~Leb(OG,/ Oc) } A(6,3) = (1/D,)[(p/7)B~(OR,/O0)
(B,/7)(ar,/ao)]
A (6,4) = (1/D~) {(p/7)[B2Le+Ba(aR,/Oc)] -- (B{7) (or,/Oc) - B~pLvo} A ( 6 , 6 ) = (l/D,)(--pLe/7)--3B~Leb)
(A20)
The steady state surface temperature 0~, is determined by
+ B4aAH (OG,/ Oc)]
- -
(A19)
(A13)
1. F~.ANK-KA~mNETSKn, D. A.: Diffu~sion and Heat Exchange in Chemical Kinetics, Plenum Pre~% 1969. 2. Mc GuIem, M. L. AND L~,VIDVS, L.: AIChE J. 11, 85 (1965). 3. W~I, J.: Chem. Eng. Sei. 20, 729 (1965). 4. B~aGEn, A. J..~r~D LAPIVUS, L.: AIChE J 14, 558 (1968). 5. Luss, D. AND L~E, J. C. M.: Chem. Eng. Sci. 23, 1237 (1968). 6. PADMANABHAN~ L., YANG, R. Y, K., XND LAPI:DUS, L.: Chem. Eng. Sei. 26, 1857 (1971).
374
HETEROGENEOUS COMBUSTION
7. HLAVACEK,V., KUBICEK,M., AND MAREK,M.: J. Catalysis 15, 17 (1969). 8. HLAVACEK,V., KUBICEK,M., ANDMAREK,M.: J. Catalysis 15, 31 (1969). 9. HLAVAC~K,V., KUSlCEK, M., AtCDMAUVE,M.: J. Catalysis 22, 364 (1971). 10. Bt,AVAOI~B:,V., KUBICEK,M., AND MAnEK,M.: Chem. Eng. Sci. 25, 1527 (1970). 11. Luss, D. AND AMU~'DSON,N. R.: Chem. Eng. Sci. 2~, 253 (1967). 12. DENN, M. M.: Chem. Eng. J. 41, 105 (1972).
13. Mc Gx~Avr, C. ~NI) SOLIMA~, M. A.: Chem. Eng. Set. 28, 1401 (1973). 14. ~VINEG_~RDN~a,D. K. ~N1) Se~IMVrz, R. A.: AIChE J. 1/~, 301 (1968). 15. KUnYLKO,L. A~'OEss;~-m(~H, R. H.: Fourteenth
Symposium (International) on Combustion, p. 1375, The Combustion Institute, 1973. 16. BEUSC~I, M., FIEGUTH, P.~ AND ~r E.: Chemie-Ingenieur Technik 4/~, 445 (1972). 17. WALZ,A, : Boundary Layers of Flow and Temperatare, MIT Press, 1969.
COMMENTS T. A. Brzustowski, University of Waterloo, Canada. You have shown that linear iustability may exist even if internal diffusion and combnstion within the particle are neglected. This result underlines the importaoee of the activation energy. IIave you considered the possible effect of two different activation energies, one for ti~e reaction within the particle and one for the surface reaction? It would also be interesting to consider the possibility of two activation energies for the one reaction as you have lumped it. I thbflr that your Eq. (33) could accommodate both of these situations.
Author's Reply. I have not analyzed the effect of two different activation energies, as you propose it, but this would be an interesting thing to do. Since a complete set of parameters from the experimental data was missing, the aim of the present calculation was to investigate the problem with the lowest number of parameters possible. This is also the reason why the effect of diffusion has been neglected in the nonlinear calculation, although it might alter the res~ilts considerably without changing though the fundamental character of ~he thermally driven oscillations.
R. H. Essenhiyh, Penn Stale University, USA. I am not entirely convinced that the author's predicted oscillations are those which we observed ia our studies (Ref. 15 of the paper). However, Dr. Peters' analysis provides au excellent basis for planning experiments that would support or test his explanation. Your comments on the possible significance of the following two points would be appreciated.
(1) We did find experimentally that the CO to CO2 reaction was apparently important in the oseillations, although it does uot figure explicitly in the analysis. Can it be said to appear implicitly so that this might be a distinguishing factor? (2) Likewise, our impressio~Lwas that the oscillations were not self-exciting but required an initial condition. Would this also be important for permitting a judgment on thc applicability of the theory?
Author's Reply. We ourselves were surprised to obtain immediate oscillatory instabilities with a single irreversible, first order heterogenous reaction when we used the parameters deduced from the experiments of Kurylko and Essenhigh)5 A refinement of the reaction scheme to include the homogeneous CO to C02 conversion in the SmTounding gas phase should yield some complementary information. Nevertheless, a summarization of homogeneous and the heterogeImous reaction rates, which was proposed in the discussion following Eq. (34) should be possible. This would allow to calculate the major effects by assuming an overall heterogeneous reaction rate as it was done. Concerning the second point I found that tim appearance of self~xited oseillations depended on the thermodynamic state of the system described by a certain number of parameters. These parameters K, ~, n, Le, fl, AH, a, and h, were assumed to be time-independent over the time range considered. In the e~perimental situation, these specific values of the parameters must be obtained in ,me way or another in order to make oscillations possible. This is probably which was meant by Prof. Essenhigh when he refers to the initial conditions which were required. On the other hand the stability of the system was not
THEORY OF HETEROGENEOUS COMBUSTION found to depend on initial conditions in the mathematical sense of the word.
F. C. Lockwood, Mech. Eng. Dept., Imperial College, England. The governing differential equations which are solved are of the (parabolic variety) boundary-layer kind. Since standard finite-difference techniques exist for the solution of parabolic differential equations which are capable of handling the complex boundary conditions and sufficiently fine grids, such as the method due to Patankar and Spalding, why was an approximate integral type of solution procedure adopted? The use of the approximate solution method renders difficult the distinction of solution-dependent effects from real physical ones.
Author's Reply. An integral method rather than a finite difference method was chosen to analyze this problem for two reasons: 1. A general algebraic stability condition such as Eq. (33), which permits one to estimate the influence of several characteristic parameters on
375
the behavior of the system, could not have been obtained by a purely numerical method. 2. The computational effort necessary to solve the system of Eqs. (1)-(13) is very high due to the non-linearities in the boundary conditions (8) and (9). I t is well known that finite difference methods are very sensitive to non-linearities in the boundary conditions. Nevertheless, since it was felt that the obtained results should be verified by a different numerical approach, a second-order finite difference method was used to calculate case B. The result obtained with 100 net points in both of the phases, when compared to Fig. 6 in the paper, showed that the finite difference calculation confirmed generally the appearance of oscillatory instabilities, while the amplitude and the frequency of the oscillations were different. The same result was obtained when 200 or even more net points were used in the calculation. Surprisingly enough, when only 50 net points were used, the finite d;~erence approximation was not accurate enough and the oscillations were damped or did not appear. Due to the necessary iterations the calculation took about 1000 times as much computer time as the calculation of case B with the integral method.
363
364
HETEROGENEOUS COMBUSTION
~ i AH'R(T,c)
the most genera[ case, it is assumed to diffuse through the pores into the interior of the sphere where it reacts internally with the solid combustible material. All of these processes are essentiaRy time-dependent and will be investigated on the basis of the balance equations in spherically symmetrical geometry. For analytic convenience, the density, thermal diffnsivity, conductivity, diffusion coefficient and heat of reaction are assumed to be constant. These restrictions could be relaxed if necessary. Governing Equations
= T..T. ~ Ts w a. X ~W
The non-dimenslonalized balance equations describing the system are: inside the particle ( 0 ~ O0
T~
/0~0 200\
1)
1
Oc /02c 2 0 c \ ~= ubl~+; g)
I1~ Z LU
(j
(1)
AH.r(O,c)
I I
es Cs + ~s
P(0, c)
;
in the surrounding atmosphere
ro
RADIUS [cm] I~G. 1. Schematic of the process of non-stationary heat and mass transfer in the combustion of a spherical particle. the application of arbitrarily chosen trim functions restricts the number of possible pertubations causing instability. Therefore the aim of the present theory is to detect causes of combustion instabilities, rather than to assure the stability of a state of a system.
(l~v< ~ )
00 --~020 200 0r Or?~ ~ 07
0
0
r
0c
(3)
_ /02c. 2 0 c \
These equations are subject to initial conditions 0(7, 0) = e0(,)
(5)
c(~, 0)= ~(n)
(6)
and to boundary conditions at the center of the sphere, the surface and at infinity: ~?=0:
The Physical Model A single spherical particle of combustible material is considered in an infinite stagnant oxidizing atmosphere (see Fig. 1). A heterogeneous reaction X I + solid--*X, that converts the oxidizing gaseous component Xt into the gaseous product X~ takes place at the surface and in the interior of the particle. The heat of reaction generated on the surface is transferred by heat conduction and radiation to the surroundings and into the interior of the spheres. The oxidizing component is supplied by diffusion from the surroundings. In
oola~=o,
(7)
octan=o
'r Q = -- k,(00/0,7) = -- (00/0'~)+ QR(0)-'2d-/" R(0, c) (S) ./= --pLcb(0c/0,7) = -- Le(0c/e,1)+ R(0, c) (9) y----~ o0 :
O= 1,
c= 1
(10)
THEORY OF HETEROGENEOUS COMBUSTION Equations (8) and (9) represent heterogeneous balance equations that relate the ~wo phases. The radiation at the surface QR(0), the reaction rate inside the particle 17(8, c) and the heterogeneous reaction rate R(O, c) at the ~rface are defined as
F = kb exp(--3b/~)'C ~
(12)
R = k exp(--fl/O~,), c. ~
(13)
where reverse reactions have been neglected. Consider now a finite disturbance 0(~, r)=O(V, r)--O~(~)
(14)
5(~, r ) = c(~, 7)-c,(~)
(1~)
of the steady state solution 0,(~/), c,(v) of the system. When they arc introduced into Eqs. (1)(10), due to their linearity Eqs. (3)-(7) and (10) remain unchanged, while Eqs. (t), (2), (g), and (9) take the fona~
a~ de _
/a~. 20~. ~H.(r-r~) /o~5 2 oe )
were chosen taking into account the fact t h a t the finite disturbance must decay at large distgnces fronl the sphere:
[~O')+dt(r)(n-- 1)
O(n, r
(11)
Q~=r
(r-r~)
+&(r)(n-1)~] 9 expl--'Kn-1)]
(18)
-re~(z)(V--1)~] 9 e x p [ - - 7 ( ~ - - l ) ]
(19)
5(~, r
The exponential decay parameter 3"= 10.O was chosen such t h a t good agreement with the linear stability conditions of Winegardner and Sclunitz ~4 was obtained in the case of their particular example. When the coefficicnts in Eqs. (16) (19) are eliminated (~e Appendix), integration of the finite disterbance Eqs. (3a) and (4a) in the form
/,
o~d,=j,f-/o~ 2o~d
(2o
/ O'~C 2~5\
(2,
-
~
(la)
(2a)
Q = - ~,( o~/ o, ) = - (o~/o~) +(QR--Q~,)--AH(R--R,)
365
(8a)
yields the integral relations that determine the solutions. The weighting function ~ was introduced here to obtain analytical integrals. When this is done inside the particle the values close to = 1 are emphasized whereas the values at, ~/---~O are neglected. For this reason the differential equations at n = 0 together with the condition of symmetry in the limit */--~0
J = -- pLeo( O~/On) = -- Le( O~/&l)+ (R - - R~)
OO/O'r= 3a(O2~/ar
(22)
(9a)
Os/a~= 3LeKaW a~~)
(23)
Since integra2 relations are to be used to study the behavior of the above system, appropriate trial functions are required to approximate the solution in the domain of integration. In the present case it was found t h a t a Pohlhausen-type polynomial trial function could approximate the finite disturbance profiles inside the particle:
were used tot the integral relations inside the particle. When this is done there results a system of two second-order equations describing the time dependency of the surface temperature and concentration, and two first-order equations for heat transfer and diffusion flux at the surface
0(~, r ) = ~.(r)q-aKr) (1-- r/)+o~(v) (l--v) ~
BI(O2~da~)+ B2(Og~/O~)+ ~(a/o~)E(Q~ - Q..)
+aK'r)(1--n)~+a~(r)(1--~) '
(16)
e(n, r) = 5~(r)-4-bKr) ( 1 - ~/)+b~(r)(1--n) 2 -I-bdr)(1--n)~+b~(r)(1--n) ~ (17) In the surrounding gas phase the trial functions
- ~.
( R - R.)~-- B~(aQ/O~-) = Q - ~
-E(Qn-QRs)-Att.(R-R.)]
(24)
Bt(O~5~/Or:)+LeB~(c3~/Or)-l-Ba(O/Or)['R--R,] --Ba(OJ/OT)=Le(J'--Le'5.,--(R-- R.)) (25)
HETEROGENEOUS COMBUSTION
366
(7a/k,) (O0./Or)-'}-O~/ Oa--}-20a( O'O~/Ov)
(aH/t,c,)(o/o~.)(r- r,),=~
-
=
60(a~/l~)Q"F60(aAH/pc~,)G
(26)
(7/0) (OJ/O'r)-I-O~6~/0~"}- 20Loo(66~/Or) + (l/p) (0/0T) (F-- F~) ~-, =
60(Le~/o).]-6O(Leb/p)G
thus stability, while a positive real part of any one of the elgcnvalucs, signifies exponential growth of the disturbance and thus instability. The latter is the ease, as may readily be seen, if the coefficient p0 in Eq. (32) takes a negative value. From this condition, using the values of the matrix A, an algebraic criterion for instability may be derived after some elementary manipulations. The system is unstable if
(27)
AH. (OR,/OO+OGJOO)(,gR,/Oc+OGg'Oc) where
< {-- 1--OQ~cJOO+AH(ORjO0 G= K(r-
r~),=~+ ~ ( r - r , ) , ~
(28)
and B~= ( ~ + 3)/~ ~
B~= (3~+6)/v 2 B~= (3~-1-6~--6)/'y4 (29)
The second-order derivatives in Eqs. (26) and (27) may be eliminated by using Eqs. (24) and (25). When each of the second order equations is replaced by two first-order ones, tlmre results a system of six first order nonlinear differential equations t h a t are to be solved simultanconsly. The unknowns of this system may be written as a vector
9 = I0., ~, oOJo~, o~Io~, Q, J}
(3o)
which represents the disturbances of the steady state solution.
Linear Stability Analysis The stability of a steady state solution for infinitely small distances may be analyzed by linearizing the nonlinear terms in Eqs. (24)-(27) around a vanishing vector of disturbances. One obtains a sixth order matrix equation of the form
Ox/Or= Ax
(31)
The elemenk~ of the matrix A (see Appendix) are evaluated for the steady state. The eigenvalues of the system (31) may be obtained by solving its characteristic equation
hn-}-p~5"~-p4M-'~-psh3+p~?--F-pl;k.4-po=O (32) Negative real parts of the eigenvalues signify damping of an infinitely small disturbance and
~-OGdO0) }. { Le-}-ORJOc-I-OG,Oc}
(33)
It is interesting to note that neither the diffusion coefficient nor the thermal diffusivity or thermal conductivity inside the particle appears in Eq. (33). t51rthennore, by combining the reaction rates in the interior of the sphere with the heterogeneous reaction at the surface, an overJdl heterogeneous reaction rate may be introduced
Re*= R,+G,
(34)
If R, is replaced by R,* with G,= 0 in Eq. (33), the stability condition takes the same form as if only a heterogeneous reaction was considered. Thus the stability condition (33) remains unchanged if a combnstion model is adopted that neglects internal diffusion and assumes the interior reaction to takc place at the surface of the sphere by summarizing the reaction rates appropriately. Note t h a t Eq. (33) is a sufficient but not a necessary condition for instability to small disturbances. Using the Hurwitz stability criterion a complete ~ t of six inequalities may be derived. Numerical investigations show however that for a large range of parameters of interest, Eq. (33) is decisive for the change from stability to instability. As it is an algebraic relation which is easily applied to a given problem, it may be valuable for practical purposes.
Nonlinear Stability Analysis The above theory has been applied to the problem of combustion instabilities reported by Kurylko and Essenhigh. 15Iu their experiments on the combustion of a radiation-ignited spherical carbon particle these authors found oscillatory instabilities with temperature excursions between 25~ and 100~ to occur at furnace tempcratures of about 1000~ In view of the above results
THEORY OF HETEROGENEOUS COMBUSTION I/.00
i
r
t i iI
~L 1300
L/J
ik/ O_
367
1200
LU
i
]
I i '~
"
J
'i /
1100
','4-I
~A
~,o
1000 0,5
1.0
1.$
2.0
RADIUS [cm] Fro. 2. Steady state solutions. it was considered t h a t diffusion into the interior and internal reaction could be neglected in the analysis of this problem. Thus J and F are zero and Eq. (27) m a y be dropped in the system (24)(27). A system of five nonlinear first order differential equations was integrated numerically b y a Runge-Kutta procedure. Consistent initial values were obtained by determining the eigens vectors of the iinearized system (31). From the experimental data given~5 the following parameters could be established and used to calculate the steady state solutions (See Appendix) : r0 = 0.635 cm
T~ = 1000~
0~ = 1.22
Le = 1.5
f~ = 24.5
~,H= 10.0
The reaction order was assumed to be n ~ 1. A major difficulty consisted in choosing the preexponential factor ]c and the radiation coefficient (r. Figure 2 shows a plot of the steady state surface temperature over the radius r0 for three different values of k. The condition 0~= 1.22 at r0= 0.635 relates a value of g definitely to each value of ]z (see Table I). In Fig. 2 linearly unstable steady states are represented by a dashed
line. The largest eigenvalues for the three cases at r0= 0, 0,,= 1.22 are also given in Table I. I t is seen from Fig. 2 t h a t in case A, three steady states exist at ru=0.635, where only the lower one is linearly stable. In case B, only one steady state exists for all values of r0 which is linearly unstable at to= 0.635, while in case C, the singlc steady state at this radius is stable. A plot of steady state heat generation by reaction and heat transfer by conduction and radiation above the surface temperature is given in Fig. 3 for case A. The steady state solutions are those where the two liues cross e~ch other. By considering the slope of both lines, it is generally argued that the upper and the lower steady state is stable while the intermediate is unstable. Since linear stability analysis predicts instability TABLE I Case
k
~
;~(max)
A B C
3.4.10 ~ 1.0-10 TM 1.0.10 ~I
9.834 11.26 13.65
270.6 15.55-{-41.96 i -8.223
Evaluated at r0 ~ 0.635 cm, 9~ = 1.22
HETEROGENEOUS COMBUSTION
368 is
//
IO n., IM u, g,) z
r < uJ -t,..
~
O
R
§ QC
CASE A
/
o
-$
J
toflo
000
WALL
1200 TEMPERATURE
I1.00
[K]
FzG. 3, Steady state hea~ balsalce in case A.
1300-
ILl n,"
120@
ILl UJ I,--
1100-
100@
TIME
'r =
tag ro2
[_]
Fzo. 4. Transient behavior of the surface temperature s t ~ t i n g from the upper steady state sohztion in case A.
T H E O R Y OF H E T E R O G E N E O U S C O M B U B T I O N
369
y
IS
10
/ 1
AH- R
~/OR*
CASE B
J
Y
-w Boo
Oc
1200
tooo
1400
WALL TEMPERATURE [K] Fx~. 5. Steady state heat balance in case B.
50
^~AAA#
,,4 Z
AI
< < -SO 5O
2,0
2.5
3.0
~.0
L.S
$.0
2 <
o. X
w
~
-SO
3.S
TIME
I"-- ~
r2
5.5
[-]
FIO. 6. Transient oscillatory behavior of the surface temperature vadation in case B.
370
HETEROGENEOUS COMBUSTION 2S
----..
r162 .p uJ
-25
/
50 uJ
~X -50 0.0,~ ,
z 0
iJ 7 o
,,.
-o,0~
3.50
3.55
3.60
TIME v = t'ro~
3.65
3.70
[-]
FIG:. 7. Surface heat transfer, temperature and concentration variations of one oscillation in case B. for the upper steady state, a nonlinear calculation was carried out to determine its transient behavior. Figure 4 shows surface temperature against time as determined by the numerical solution of Eqs. (24)-(27). As the tcmpcrature disturbance was assumed to be only 0.03~ an induction period was required before the temperature dropped rapidly. The chemicaI reaction is quenched and aftcr an overshoot %o lower surface temperatures the stable steady state at 1042~ is reached at r = 3. This calculatioa suggests that the experimentally observed oscillations do not correspond to a situation where 3 steady states exist, since in this case the lower stabIe state would be reached. Furthermore it is shown that physical argu-
mentation may fail in certain cases while linearized mathematical analysis may predict tile correct physical behavior. If the steady state heat balance is considered for case B (Fig. 5), the slopes of the two curves differ considerably at the siugle steady state. Nevertheless the nonlinear calculations (Fig. 6) show exponentially growing oscillations after an induction period of r = 1.5 at ~he assumed initial temperatm'e pertnbation at 0.014~ At r = 4.0 the oscillations reach a constant limiting cycle with an upper amplitude of 52~ and a lower limit of 40~ To aoalyzc the oscillations, the cycle between r = 3.5 and r = 3.7 has been considcred in detail. Figure 7 shows the time-dependent behavior of heat trausfer, surface tem-
THEORY OF HETEROGENEOUS COMBUSTION
371
SO
" ~
T ~ , 3.62 "'
~ /
0
-so
.• o
T = 3.Se /,"-"--
t
2 RADIAL DISTANCE ;7 = r/r o
v
9 3,67
-SO
RADIAL DISTANCE -- r/r o
Fro. 8. SpaciM distribution of the temperature variation of one oscillation ia case B. perature and surface concentration. There is a phase shift between heat transfer to the surface and the su~ace temperature. Heat is continued to be furnished from the iuside of the particle at a nearly constant rate between r = 3.52 and v = 3.57 while tlie temperature drops continuously. Then, when the temperature begins to rise again, the heat transfer suddenly changes its sign. Heat is conducted at a high rate into the interior of the particle and thus stops the rise of the surface temperature. Prachcally no phase shift is found between the surface temperature and conceutratiou. These obscrvations lead to the conclusion that the experimentally observed oscillations are caused by thermal instabilities ratber than by internal diffusion. The temperature profiles inside and outside the particle for the same period of time are shown in Fig. 8. The amplitude maximum lies inside the particle at a certain distance from the surface. Heat exchange by conduction between this region
and the surface causes a time lag which is responsible for the instabilities. For the above calculation, values of ]c~=2.2 and a = 1.0 were asstoned. Since the value of thermal conductivity k~ is not needed to determine the steady state, it may be altered without changing the above considerations. When values of kt = 2.4 and higher were chosen, it was fmmd that the steady state was stable to finite disturbances. This indicates that the influence of porosity on the stability limits found by Kurylko and Essenhigh~b may be indirectly due r a change in thermal conductivity. This may explain the stability at low porosities. Case C was found to be stable to finite disturbances for all values of k,. From the comparison of the values of k and a iu the three cases considered, where only the intermediate one exhibits the possibility of oscillatory instabilities, it is concluded that there is only a narrow range of
HETEROGENEOUS
372
parameters in which the observed oscillations could occur. This is in agreement with the observations of Kurylko and Esseahigh ~ and may explain those limitations to the region of instability which were not investigated in the present analysis.
COMBUSTION
B r 7 9
6 Conclusions
1. Linearly unstable sLeady states may exist in the heterogeneous combustion on a spherical particle even if internal diffusion and combustion are neglected. 2. To detect the eau~s of these instabilities, nonstationary heat conduction and diffusion to the surroundings, i.e., non-stationary Nusselt and Sherwood munbers, must be considered. 3. If multiple steady states exist for a given problem, the assmnptiou t h a t the outer ones are always stable may not be valid. 4. Experimentally observed oscillatory ins~bilities can be explained b y a phase shift between the surface temperature and heat transfer from the interior. Nomenclature
a ab/% dimensionless thermal diffnsivity c cl/c1| mass fraction of component 1 cp c#/%g dimensionless specific heat E G
activation energy combined internal reaction rate Eq. (28) 511 (hl--l~)ct~/cm T~ dimensionless beat of reaction h enthalpy J j~ro/p#aocl~dimensionless diffusion flux k koaro(p~l~)"-l/aa dimensionless preexponential factor of the heterogeneous reaction kb "]~r dimensionless preexpcnential factor of the internal reaction kt ktb/kta dimensionless thermal conductivity Le D./a# gas phase Lewis number Le~ D~/% internal diffusion Lewis number n reaction order Q q~ro/k,~ T~ dimensionless heat transfer to the surface Qc 0 ~ - 1 dimensionless convective heat transfer Q~ dimensionless heat transfer by radiation Eq. (11) R a universal gas constant R dimensionless heterogeneous reaction rate Eq. (13) R* dimensionless overall heterogeneous reaction rate Eq. (34) r radial distance r0 radius of the sphere T Temperature t time
X
P o~t
r
E/RgT~ dimensionless activation energy dimensionless internal reaction rate Eq. (12) decay parameter emissivity for radiative heat loss r/ro dimensionless radial distance TIT** dimensionless temperature eigenvalue Pb/Pgdimensionless density ~teroT~/k,adimensionless radiation coefficient Boltzmann radiation constant agt/ro~ dimensionless time
Supexscript disturbance from the steady state
-
Subscript b q s w ~r
solid body surrounding ge~ steady state surface (wall) at infinity Appendix
The coefficients in Eqs. (16)-(19) are eliminated by using the boundary conditions at the origin and at the interface as well as by the requirement t h a t the trial functions mus~ fulfill the differential equations at 7 = 1. Finally, Eqs. (16) and (17) must be symmetrical around the origin, One obtains the following equations: al+2~+3az+4ad= 0
(All
l)1+2b2+ 3b3+ 4b4= 9
(A2)
~ = Q/k,
(A3)
b~= J/pLeb
(A4)
o~/o~=a(2ar 2a,)+(~/z~p)(r- r.)~.l (AS) o~/or= L~(2b2--2h)--(1/p)(r--r.)~.l (A6) a~-{-2a2-1-4 a s + 8a4= 0
h+~b~+4bs+Sbd= O
O=--(d,-7~)+(Q~--Qn~
(A7)
(AS)
R.) (AS)
] = -- Le(et-Te.)+ (R-- R,)
(A10)
O~/Or=26z+2d~(1--7)+~(7~--27)
(All)
0~/0~= Le[2~+2~(1--v)+~(v'--2x)] (A12)
T H E O R Y OF H E T E R O G E N E O U S COMBUSTION The elements of the matrix A, Eq. (31) have the following form: A(1,1)= A ( 1 , 2 ) = A ( 1 , 4 ) = A(1,5)= A(1,6)= 0.0
373
where
Bi = 60B{7,
B5 = 20B3/7,
B 6 = 20B1/7
A (2,1) = A (2,2) ~ A (2,3) ~ A (2,5) = A (2,6) = 0.0 A (3,6) = A (4,5) = A (5,6) = A (6,5)= 0.0
(A14) and
A(1,3) = A (2,4)= 1.0
D,=BI+k~Bj7a,
A (3,1)= ( l / D , ) { - - 1-- (0/oM)
D,=B1+pB~/7
(A15)
9[QR,--AH(R,+B4aG,)]}
The internal reaction rate F, in Eq. (A13) is evaluated at y = 1. The steady state solution for n = 1, F = 0 a n d J = 0 is given b y
A ( 3 , 2 ) = (1/D,) { AII (a/Oe)[ R,+ B4aG,]} A (3,3) = (1/Dr) { -- B2-- s3(a/ao)
9[QR,--AH(R,+F,/7)]--B~k,} A(3,4) =
O/D,)[Ba&H(a/ac)( R , + r J 7 ) ]
c,(~)=BJ~+l
A(3,5) = (1/Dt)(1--B4a)
(1_>,/>oo)
(A16)
where
a ( 4 , 1 ) = (1/D~) {-- (a/ao)
9(LeR~+B~LebG~) }
B~:
k- exp(-~/o~,) Le+k. exp(--fl/8~,)
(A17)
A ( 4 , 2 ) = (1/De)I-- Le~- (O/Oc)
O,(*/)= 0~,
9(LeR~+B4LebG~)] A (4,3) = (1/D~)[-- B3(O/00) (R,+ rdT)]
(0_>n_> 1)
(A18)
and
A (4,4) = (I/D,) [-- B2Le-- Ba(O/Oc) 0,=Bdn+l
9(R,+ r d 7 ) - B~pLeb] A ( 4 , 6 ) - (1/D~) ( Le-- B4Leb) A (5,1) = (1/D,)I (k,/7a)[l+ (a/ao)
( l > y > ~o)
where Bt= 0~-- 1
9 ( Q n , - &H R,)_-]+ B#AH(aG,/00) } A (5,2) = (1/Dr) I - - (kJ7a)AH(OR,/Oe) A(5,3) = (1/D~) {(k,/7a)[B~+B,(a/ao) x ( Q R , - A H R , ) ] + (B1/7)AH(OF,/OO)- Bdv,}
0~,-- l + a ( 0 ~ , ' - - 1)
AHLek exp(--fl/8~) ~
A(5,4) = ( l / D , ) [ - - (kd7a)BaAH(OR,/ae)
= o
(A21)
+ (B~/7)&H(aFo/Oc)] A (5,5) = (lIDs) (-- k~/Ta-- 3B~a)
REFERENCES
A(6,1) = (1/D~)[(pLe/7)(OR,/O0)
-- B4Leb(OG,/aO)] A(6,2) = (1/D~)I (p/7)[Le~+ Le(OR,/Oc)]
-- B~Leb(OG,/ Oc) } A(6,3) = (1/D,)[(p/7)B~(OR,/O0)
(B,/7)(ar,/ao)]
A (6,4) = (1/D~) {(p/7)[B2Le+Ba(aR,/Oc)] -- (B{7) (or,/Oc) - B~pLvo} A ( 6 , 6 ) = (l/D,)(--pLe/7)--3B~Leb)
(A20)
The steady state surface temperature 0~, is determined by
+ B4aAH (OG,/ Oc)]
- -
(A19)
(A13)
1. F~.ANK-KA~mNETSKn, D. A.: Diffu~sion and Heat Exchange in Chemical Kinetics, Plenum Pre~% 1969. 2. Mc GuIem, M. L. AND L~,VIDVS, L.: AIChE J. 11, 85 (1965). 3. W~I, J.: Chem. Eng. Sei. 20, 729 (1965). 4. B~aGEn, A. J..~r~D LAPIVUS, L.: AIChE J 14, 558 (1968). 5. Luss, D. AND L~E, J. C. M.: Chem. Eng. Sci. 23, 1237 (1968). 6. PADMANABHAN~ L., YANG, R. Y, K., XND LAPI:DUS, L.: Chem. Eng. Sei. 26, 1857 (1971).
374
HETEROGENEOUS COMBUSTION
7. HLAVACEK,V., KUBICEK,M., AND MAREK,M.: J. Catalysis 15, 17 (1969). 8. HLAVACEK,V., KUBICEK,M., ANDMAREK,M.: J. Catalysis 15, 31 (1969). 9. HLAVAC~K,V., KUSlCEK, M., AtCDMAUVE,M.: J. Catalysis 22, 364 (1971). 10. Bt,AVAOI~B:,V., KUBICEK,M., AND MAnEK,M.: Chem. Eng. Sci. 25, 1527 (1970). 11. Luss, D. AND AMU~'DSON,N. R.: Chem. Eng. Sci. 2~, 253 (1967). 12. DENN, M. M.: Chem. Eng. J. 41, 105 (1972).
13. Mc Gx~Avr, C. ~NI) SOLIMA~, M. A.: Chem. Eng. Set. 28, 1401 (1973). 14. ~VINEG_~RDN~a,D. K. ~N1) Se~IMVrz, R. A.: AIChE J. 1/~, 301 (1968). 15. KUnYLKO,L. A~'OEss;~-m(~H, R. H.: Fourteenth
Symposium (International) on Combustion, p. 1375, The Combustion Institute, 1973. 16. BEUSC~I, M., FIEGUTH, P.~ AND ~r E.: Chemie-Ingenieur Technik 4/~, 445 (1972). 17. WALZ,A, : Boundary Layers of Flow and Temperatare, MIT Press, 1969.
COMMENTS T. A. Brzustowski, University of Waterloo, Canada. You have shown that linear iustability may exist even if internal diffusion and combnstion within the particle are neglected. This result underlines the importaoee of the activation energy. IIave you considered the possible effect of two different activation energies, one for ti~e reaction within the particle and one for the surface reaction? It would also be interesting to consider the possibility of two activation energies for the one reaction as you have lumped it. I thbflr that your Eq. (33) could accommodate both of these situations.
Author's Reply. I have not analyzed the effect of two different activation energies, as you propose it, but this would be an interesting thing to do. Since a complete set of parameters from the experimental data was missing, the aim of the present calculation was to investigate the problem with the lowest number of parameters possible. This is also the reason why the effect of diffusion has been neglected in the nonlinear calculation, although it might alter the res~ilts considerably without changing though the fundamental character of ~he thermally driven oscillations.
R. H. Essenhiyh, Penn Stale University, USA. I am not entirely convinced that the author's predicted oscillations are those which we observed ia our studies (Ref. 15 of the paper). However, Dr. Peters' analysis provides au excellent basis for planning experiments that would support or test his explanation. Your comments on the possible significance of the following two points would be appreciated.
(1) We did find experimentally that the CO to CO2 reaction was apparently important in the oseillations, although it does uot figure explicitly in the analysis. Can it be said to appear implicitly so that this might be a distinguishing factor? (2) Likewise, our impressio~Lwas that the oscillations were not self-exciting but required an initial condition. Would this also be important for permitting a judgment on thc applicability of the theory?
Author's Reply. We ourselves were surprised to obtain immediate oscillatory instabilities with a single irreversible, first order heterogenous reaction when we used the parameters deduced from the experiments of Kurylko and Essenhigh)5 A refinement of the reaction scheme to include the homogeneous CO to C02 conversion in the SmTounding gas phase should yield some complementary information. Nevertheless, a summarization of homogeneous and the heterogeImous reaction rates, which was proposed in the discussion following Eq. (34) should be possible. This would allow to calculate the major effects by assuming an overall heterogeneous reaction rate as it was done. Concerning the second point I found that tim appearance of self~xited oseillations depended on the thermodynamic state of the system described by a certain number of parameters. These parameters K, ~, n, Le, fl, AH, a, and h, were assumed to be time-independent over the time range considered. In the e~perimental situation, these specific values of the parameters must be obtained in ,me way or another in order to make oscillations possible. This is probably which was meant by Prof. Essenhigh when he refers to the initial conditions which were required. On the other hand the stability of the system was not
THEORY OF HETEROGENEOUS COMBUSTION found to depend on initial conditions in the mathematical sense of the word.
F. C. Lockwood, Mech. Eng. Dept., Imperial College, England. The governing differential equations which are solved are of the (parabolic variety) boundary-layer kind. Since standard finite-difference techniques exist for the solution of parabolic differential equations which are capable of handling the complex boundary conditions and sufficiently fine grids, such as the method due to Patankar and Spalding, why was an approximate integral type of solution procedure adopted? The use of the approximate solution method renders difficult the distinction of solution-dependent effects from real physical ones.
Author's Reply. An integral method rather than a finite difference method was chosen to analyze this problem for two reasons: 1. A general algebraic stability condition such as Eq. (33), which permits one to estimate the influence of several characteristic parameters on
375
the behavior of the system, could not have been obtained by a purely numerical method. 2. The computational effort necessary to solve the system of Eqs. (1)-(13) is very high due to the non-linearities in the boundary conditions (8) and (9). I t is well known that finite difference methods are very sensitive to non-linearities in the boundary conditions. Nevertheless, since it was felt that the obtained results should be verified by a different numerical approach, a second-order finite difference method was used to calculate case B. The result obtained with 100 net points in both of the phases, when compared to Fig. 6 in the paper, showed that the finite difference calculation confirmed generally the appearance of oscillatory instabilities, while the amplitude and the frequency of the oscillations were different. The same result was obtained when 200 or even more net points were used in the calculation. Surprisingly enough, when only 50 net points were used, the finite d;~erence approximation was not accurate enough and the oscillations were damped or did not appear. Due to the necessary iterations the calculation took about 1000 times as much computer time as the calculation of case B with the integral method.