Simultaneous higher order reactions and diffusion in a slab


Material Information

Simultaneous higher order reactions and diffusion in a slab
Physical Description:
x, 211 leaves : ill. ; 29 cm.
Koopman, David, 1954-
Publication Date:


bibliography   ( marcgt )
theses   ( marcgt )
non-fiction   ( marcgt )


Thesis (Ph. D.)--University of Florida, 1992.
Includes bibliographical references (leaves 208-210).
Statement of Responsibility:
by David Koopman.
General Note:
General Note:

Record Information

Source Institution:
University of Florida
Rights Management:
All applicable rights reserved by the source institution and holding location.
Resource Identifier:
aleph - 001862353
notis - AJT6820
oclc - 28922401
System ID:

Full Text







This dissertation is dedicated to my parents and sister without whose support and

encouragement it would never have been possible.


I would like to gratefully acknowledge the assistance and inspiration provided by

Dr. Hong H. Lee both through his constant assistance in this endeavor and through his

textbook, Heterogeneous Reactor Design. I would also like to express my appreciation

to Rutherford Aris, who collected and reviewed many of the early accomplishments in

this area of research. Without the accomplishments of those who have gone before, this

work would not have been possible.


ACKNOWLEDGEMENTS .................................. iii

KEY TO SYMBOLS ..................................... vi

ABSTRACT .......................................... ix

Introduction ....................................... 1
General Problem Statement ............................. 4
The Zeroth-Order Reaction ............................. 8
The First-Order Problem .............................. 10
Effectiveness Factors ................................. 17
The Generalized Thiele Modulus ......................... 22
Shape Factor Normalizations ............................ 26
Summary and Goals ................................. 27

Overview of Second-Order Reaction and Diffusion .............. 30
Governing Equations and Dimensionless Groups ................ 34
Solution of the Dirichlet Problem ....................... 45
The Two Fundamental Parameters ..................... 52
The Weierstrass Pe-Function, .p .......................... 66
Concentration Profiles From the Analytical Solution ............. 75
First-Order Subcases of the General Problem .................. 78
Internal Effectiveness Factors of Second-Order Reactions .......... 84
Concluding Remarks ................................. 104

Introduction to Third-Order Reaction and Diffusion ............. 110
Three-Parameter Equations and Concentration Profiles ............ 111
Restrictions on the Feasible Parameter Space .................. 126

Internal Effectiveness Factors of Third-Order Reactions ........... 132
Concluding Remarks ................................ ..... 144

Introduction ........................................ 146
The Robin's Problem ................................ 148
The Single Species Robin's Problem ....................... 155
Autocatalytic Reactions .............................. 164
Multiple Reaction Systems ............................. 171
Concluding Remarks .................................. 182

Recommendations .....................
Future W ork .........................


. 185

. 189
. .. .. 189

REFERENCE LIST ................... ................... 208

BIOGRAPHICAL SKETCH ................................ 211



a coefficients of the single species second-order kinetic expressions

A coefficients of the Taylor series expansion of q

Bi Biot number for mass transfer

c' normalized dimensionless concentration

c" normalized dimensionless concentration

c dimensional concentration

C expansion coefficients of the Weierstrass elliptic Pe-function

D effective diffusivity

g polynomial invariants for the Weierstrass Pe-function

k reaction rate constants

L characteristic length of the slab (i.e. half width)

n order of a power law rate expression

p coordinate system identifier for slab, cylinder or sphere

p Weierstrass elliptic Pe-function

q derivative of a concentration with respect to a length scale variable

r kinetic rate expression

x dimensional position variable

z dimensionless position variable, x/L

Greek Letters

a constant in the integration of q, also dummy integration variable

S arbitrary constant in the development of Z(0)/c,

v stoichiometric coefficient

t Thiele, or first-order, or generalized, modulus

?7 effectiveness factor

I' second-order modulus

a transformation variable during integration

7 transformation variable during integration

X third-order modulus


A specific chemical species

B specific chemical species

C specific chemical species

eq equilibrium value

f forward direction

G generalized, or normalized, modulus

i property appropriate to generic chemical species, I

int internal

j,k properties appropriate to generic chemical species, J and K

l,m properties appropriate to generic chemical species, L and M

m+ associated with the surface of repeated equilibrium roots based on the positive

square root

m- associated with the surface of repeated equilibrium roots based on the negative

square root

o at the midplane, x = 0 or z = 0

P specific chemical species; or related to product formation in the rate r

Q specific chemical species

r reverse direction

R specific chemical species

s at the surface of the slab medium, x = L or z = 1

1st first-order

Abstract of Dissertation Presented to the Graduate School
of the University of Florida in Partial Fulfillment of the
Requirements for the Degree of Doctor of Philosophy



David Koopman

December 1992

Chairperson: Dr. Hong H. Lee
Major Department: Chemical Engineering

Analytical solutions are presented for a broad class of simultaneous steady-state,

isothermal second-order and third-order reactions and diffusion in a slab-like medium

with constant diffusivities. All intrinsically third-order reversible and irreversible kinetic

schemes can be unified into a single three-parameter problem, while similar second-order

kinetic schemes can be unified into a single two-parameter problem. Analytical solutions

for dimensionless concentration profiles, internal effectiveness factors, etc. are

determined as functions of these dimensionless parameters: a first-order, or Thiele,

modulus, ?,,, a second-order modulus, I, and a third-order modulus, X, for third-order

reactions only. The analytical solutions employ the Weierstrass elliptic Pe-function, p,

to characterize the spatial dependence of the concentration profile. The analysis provides

semi-quantitative conditions for both diffusion free and diffusion limited behavior. The

modulus, 9, is especially adjustable by varying physical parameters. Other limiting

interrelationships between the three moduli are developed that restrict the size of the

feasible solution space for conventional second-order and third-order reactions. Existing

approximation methods are only partially validated. Errors exceeding 20% are quite

possible. Improvements to some methods are suggested by the analysis. Autocatalytic

second-order and third-order reactions do not conform to established norms for

conventional reversible and irreversible reactions. Internal effectiveness factors well in

excess of unity are observed. The analytical solution methods developed here are valid

for either the Dirichlet or Robin's boundary condition at the surface. These methods

apply to varying degrees for selected multiple reaction systems, such as parallel reactions

of different orders and Van der Vusse kinetics. As an overall consequence of these

developments, the number of elementary chemical reactions possessing an analytical

solution for simultaneous reaction and diffusion has increased tremendously.



The study of the theory of isothermal, coupled, steady-state reaction and diffusion

problems dates back to the time of E. W. Thiele and his seminal paper: "Relation

Between Catalytic Activity and Size of Particle". The concepts of a modulus,

characteristic of the reacting species, and of a ratio of actual achieved rate to maximum

possible external rate of a chemical reaction were presented there. Such moduli bear

Thiele's name today, while the reaction rate ratios are referred to as internal effectiveness

factors. Thiele derived solutions for the case of a first-order irreversible single species

reaction with constant diffusivity and a Dirichlet surface boundary condition for both a

sphere and an infinite slab of finite thickness. Also derived was an expression for the

internal effectiveness factor of a single species irreversible second-order reaction in a

slab. More will be said about this result later.

The fifty years following Thiele's work saw the addition of analytical solutions

for single reactions with zeroth-order kinetics, with Robin's boundary conditions, in

cylindrical coordinates, and a thorough analysis of the first-order problem with multiple

reactions, as well as geometries with two or three characteristic length scales such as the

rectangular parallelepiped. It was recognized by Aris (1957, 1975) that use of the proper

characteristic length scale brought the results for different geometries into asymptotic

agreement. It was further recognized simultaneously by Bischoff (1965), Aris (1965),

and Petersen (1965) that the use of a specially normalized, or generalized, Thiele

modulus brought the internal effectiveness factor results for arbitrary kinetics into

asymptotic agreement at both large and small modulus magnitudes. This work was

collected and summarized by Rutherford Aris in 1975. Some of the emphasis in this

field shifted from analytical to numerical solutions with the proliferation of high speed

computers in the 1960's. The first studies on the behavior of Langmuir-Hinshelwood

kinetics coupled with diffusion began to appear; see Aris (1975) or Satterfield (1970) for

summaries. Problems with complex pore networks, nonisothermal reactions, variable

diffusivities, reactions with volume changes, etc. were analyzed by numerical methods.

However, after all of this effort, published analytical solutions for the

concentration profiles of single reaction systems, either reversible or irreversible, existed

only for the following rate laws:

(1) r, = k,, c, > 0 rf = 0, c, = 0

(2) rf = kc1

(3) rf = k1c, k2 restrictions on k2 may apply

(4) rf = kc, k2C2

(5) r, = kc6

where the result published for (5) was only for the concentration at the slab midplane.

This work derives analytical solutions for 25 elementary kinetic rate laws for irreversible

and reversible third-order and second-order reaction systems. It also explains how to

derive analytical solutions for five second-order and third-order reactions reversible with

a zeroth-order reaction. Additional solutions for several simple parallel reactions with

second-order and/or third-order terms are also derived, as well as solutions for about 17

types of autocatalytic behavior undergoing net second-order or third-order kinetics.

These analytical solutions are restricted to assumptions of constant diffusivity, infinite

slabs, and isothermal behavior; however, both the Dirichlet and Robin's problems are

successfully analyzed.

Past studies of problems of simultaneous reaction and diffusion are now a part of

traditional chemical engineering. Advances in this area tend to inherit some of this

established relevance. What other reasons are there to study such problems? First and

foremost, simultaneous reaction and diffusion problems have considerable industrial

significance due to the importance of supported catalysts in many commercial processes.

Secondly, while this problem was first identified over 50 years ago, there are only a

handful of single reaction kinetics which have an analytical solution in some geometry.

Methods were developed for approximating behavior in arbitrary catalyst geometries once

a solution was found in one geometry. This was an inducement to create a library of

analytical solutions in some geometry that has never been fully exploited. Thirdly,

certain types of theoretical analyses require an analytical solution as a starting point. For

example, when the reaction occurring is exothermic, there is a potential for the reaction

and diffusion problem to possess multiple steady states. One type of analysis which

assumes an isothermal pellet with a temperature gradient confined to an external

boundary layer is often valid. The analysis of the multiplicity of steady states as a

function of the problem data has been done for such first-order reaction problems,

Burghardt and Berezowski (1990). It is now possible to extend these analyses to second-

order and third-order problems. Finally, the behavior of many problems can be

approximately preserved through linearizations of the variables. The analysis of the

linearized problem is generally easier to perform. Nonlinear behavior is, nevertheless,

best studied directly. Elementary reactions involving autocatalysis are an excellent

example. They will be seen to exhibit behavior well outside the range found in linear

problems of simultaneous reaction and diffusion.

General Problem Statement

The governing equation for one-dimensional, isothermal, steady-state,

simultaneous reaction and diffusion in an infinite slab of finite thickness for the

concentration of species, i, ci, takes the form

D = vr, 1.1

where rp is the rate of product formation, Di is the effective diffusivity of species i, x is

the length scale variable measured from the slab midplane to the surface, and vi is the

stoichiometric coefficient of i in the reaction described by rp. Two geometries for the

problem in a slab are illustrated in Figure 1.1, parts a and b. Identical external

conditions are assumed when there are two exposed surfaces in order to give a symmetric





0 0
/ 0 o10

0 \ / / 0

S 0 0
0 0 v 0
0 0 0

4 4



Open Cylindrical Pores

S Sealed Walls--

x=-L x=O x=L


Figure 1.1 The Geometries for the Simultaneous Reaction and Diffusion Problem. (a)
the sealed slab, (b) the symmetric slab, (c) the sealed cylindrical pore, and (d) the
symmetric cylindrical pore.

.. Slab
w- X


problem about the slab midplane. The slab medium is infinite in the y and z directions.

Two additional geometries, or interpretations, are shown in Figure 1.1, parts c and d.

The cylindrical pore interpretation requires an assumption to neglect radial concentration

gradients, which, for surface driven heterogeneous reactions amounts to making the

kinetics homogeneous within the pore as well. The analysis to follow for second-order

and third-order kinetics requires identical conditions at +L, where L is the slab half-

width. The imposition of net gradients across the 2L wide systems gives a slightly

altered final equation to integrate, but analytical solutions may still be possible; see Ince

(1927). The boundary conditions for the Dirichlet problem are

ci = cs at x = L 1.2

-= 0 at x = 0. 1.3

The boundary conditions for the Robin's problem are

ki(if is) = Di at x = L 1.4

d- = 0 at x = 0 1.5

where C is the concentration of i in the free stream fluid away from the surface, and kgi

is the mass transfer coefficient for species i from the free stream fluid to the slab surface.

The following additional result is very important for systems where rp depends on the

concentration of more than one species

Did2 D d2c
-_- = r 1.6
v, dx2 vj dx2

Integration of eq. (1.6) from x = 0 to x gives

D, dc, D_ d 1.7
vi dx vj dx

using the boundary condition at x = 0 common to the Dirichlet and Robin's problems.

Integration of eq. (1.7) from x to x = L gives

Di D.
.-( c,-),- = "(-1.8
v v.

for the Dirichlet or Robin's problem. Eq. (1.8) can be used to eliminate the cj in rp in

favor of a single ci. Eqs. (1.8) and (1.4) can be combined to eliminate all but one of the

unknown surface concentrations in the Robin's problem boundary condition; however,

this is not necessarily desirable. (As an alternative the surface concentrations could be

carried along as implicit parameters or replaced with the internal effectiveness factor.)

This point will be examined further later in this chapter as well as in Chapter Four.

Relationships such as eq. (1.8) permit expressing the in terms of one c,, either

explicitly for the Dirichlet problem, or implicitly through the unknown c,'s for the

Robin's problem. The single species reaction rate expression that results is valid within


the slab only, since that is the range of eq. (1.8). The resulting rate expression does not

give the bulk fluid reaction rate correctly except for special cases.

The Zeroth-Order Reaction

A zeroth-order reaction is one that is independent of concentration so long as

reactant remains present. This gives a two part governing equation

Dd2 k > 0
dx2 0 e = 0

The solution for the Dirichlet problem is

S- kL2
c, 2 2D() L] 2Dc,


or, for kL2/2D, > 1, is

i ----k L-2 x* kLL2D

C -
-= = 0O x x *

x* = L 1 -s
1 kL

x* x s L 1.11



Solutions of this type are somewhat awkward because of the transition point, x', at which

both c and dc/dx are zero, superseding the usual boundary condition at x = 0. These

solutions do suggest likely forms for a dimensionless concentration and length, as well


as one additional dimensionless group, kL2/2Dc,, which is the square of the Thiele

modulus for a zeroth-order reaction, 4k. The above solutions could be written

c = $k(1 -2) Z i 1.14

S kL2 x c c 1.15
z -C
2Dc,' L' C,

or, for 4k > 1,

c = (1 2) -2(0- Q)(1 -z) z* _z < 1 1.16

c =0 Osz z* 1.17

Z = i ).

It should be noted that the Thiele modulus contains the surface concentration. This is

typical of all but first-order reaction systems, and complicates the transition from the

Dirichlet to the Robin's boundary condition problem. The solution to that zeroth-order

problem will not be needed in what follows, but zeroth-order Dirichlet problem results

will occasionally be drawn on for comparison purposes. This review also serves to

remind the reader of the possible complications that a zeroth-order reaction could bring

to a reversible reaction problem where the opposing reaction is not zeroth-order.

Generally, a real species should be associated with the zeroth-order reaction, and

safeguards should be taken to see that the concentration of this species is never less than


The First-Order Problem

First-order kinetics must be discussed in this work, since second-order and third-

order kinetics were found to behave identically to first-order kinetics for some special

cases or under asymptotic conditions. Somewhat more involved expressions for the

dimensionless concentration will be developed for second-order and third-order kinetics

than are normally used for the first-order problem. The first-order kinetics problem

review will be made using a new dimensionless concentration that was found to be the

most generally useful

c' (z) = 1.19
viL2 rp

where rp, is the rate of product formation evaluated at the surface of the slab, and z is

x/L. Note that if the reaction is in a state of equilibrium at the surface, this equation

would call for division by zero. This dimensionless concentration in eq. (1.19) is a

universal, or species independent, dimensionless concentration, since

(i ci Di ( j) 1.20
Vi v.

for the symmetric Dirichlet and Robin's problems discussed here. Using this

dimensionless concentration here will familiarize the reader with some of its attributes.


Throughout this work, reactant species will be taken as A, B,... and product species will

be taken as P, R,... The "forward" reaction rate constant, kf, will be associated with

reactants, and the "reverse" reaction rate constant, k,, will be associated with products.

Whether the forward reaction rate is actually larger than the reverse reaction rate will be

irrelevant. Consider the reversible first-order reaction A + P. The governing equations


d2k 1.21
DA d2= -krp + k 1.21

D, 2 = k, -kfA 122

Taking z = x/L along with eqs. (1.19) and (1.20) converts both eqs. (1.21) and (1.22)


d2 kfL2 krL2 c,-. 1.23
dx2 DA D

Setting either k, or kf to zero recovers the irreversible reaction problem for A or P

respectively. The minus one term in eq. (1.23) does not vanish for irreversible

reactions. The definition of c" in eq. (1.23), however, alters through rp, changing

when the problem is made irreversible. The formulation in eq. (1.23) retains the

boundary condition at x = 0 as

dc "

z =0


and has the second boundary condition

c" =0



because of the definition of c" in terms of (ci, c). Achieving a zero valued surface

concentration will give essential flexibility when trying to limit the parameter space

dimension for second-order and third-order reaction systems. Rewrite eq. (1.23) as


d2C,, 2
dc c 1.
dx = c

The square of the first-order, or Thiele, modulus, 1s,,, is recognized as the sum of the

squares of the two Thiele moduli for the irreversible forward reaction, A -, and the

irreversible reverse reaction, P -, given by


OA = L fA

The stoichiometric coefficients were assumed to be + 1, but this restriction can be lifted

in which case

OP = L .
\ p

Lt = vP ')p VA A 1.28

Expressions that possess the feature of a stoichiometric coefficient weighted sum over all

species, such as this one for the coefficient of c", will be found to persist for second-

order and third-order reaction systems. The solution of eq. (1.26) is straight forward by

conventional techniques for linear second-order differential equations with constant

coefficients and gives

cosh(cQO,,) cosh(l,,z)
c"(z) = 1.29
l acosh( ,t)

There is a less conventional technique for solving eq. (1.26), that begins by multiplying

by dc"/dz and integrating from z = 0 to z that yields

,2 2
1dc" = s(c"2 "(02) -(c" c"()) 1.30
2 dz ) 2

where c"(0) is an as yet unknown value of the solution. From eq. (1.30) one can

derive the following

dcz = 2 c"(0) 1 ,c"(0)2. 1.31
dz )z.d


The dimensionless surface concentration gradient depends only on the Thiele modulus

and dimensionless midplane concentration. Formulae of this type will recur in Chapters

Two and Three for second-order and third-order kinetics respectively. The integration

of eq. (1.30) can be continued by conventional methods after separation of variables,

ultimately leading to the expression for c"(z) obtained above.

Since c" was taken as zero at the surface, it might be thought that c"(0) is

unconstrained in magnitude. This is not the case, however, for first-order reversible and

irreversible reactions as shown in Figure 1.2. The expression for c"(z=0) can be

written as

1 1/cosh(O1st)
c "(0) = 1.32

with the two limiting values

1 1
lim c"(0) lim c"(0) = 1.33
,,,-o 2 1,,- Is

The limit at infinity is a zero of the right-hand side of eq. (1.26). It corresponds to

kcA(0) k,p(0) = 0, as can be shown by substituting into eq. (1.19). The right-hand

side of eq. (1.26) is a diffusion modified reaction rate expression, and a zero of this is

a statement of positional equilibrium in the slab medium. The limit as 4~, goes to zero

has physical significance as well. In eq. (1.26) take 4,,c" 1 = -1, the dimensionless



surface reaction rate. Integrate twice to obtain

c"(z) = + az + b.

Applying the Dirichlet boundary conditions gives

c"(z) = 1(1 2)

c"(O) = -.

Thus c"(0) = 1/2 symbolizes a uniform reaction rate from the slab surface to the

midplane. Figure 2 further suggests that 0.2 < 1,st < 5 contains the majority of the

transition from one limit to the other.

The Robin's boundary condition transforms to

Bic d" c
Biic" dc
fl dz

z= 1


where Bi, = kgL/D, is the Biot number for mass transfer for species i. Equivalently by

eq. (1.31),


Bimc",i = 2c"(O) 4c<"(0)2.

Reasons for choosing the minus sign are essentially physical. The dimensionless surface

concentration, zero, is expected to lie between the dimensionless midplane concentration,




positive, and the dimensionless free stream concentration. Therefore the dimensionless

free stream concentration is expected to be negative. A further reason for choosing the

minus sign in eq. (1.31) will be given below. Insertion of the c"f expression into eq.

(1.19) gives ci, for known cf. Knowledge of ci, plus c"(z) gives ci(x). Since ,,t is

independent of c6,, c,, etc. for first-order reactions, the final solution for the Robin's

problem is immediately obtained from the Dirichlet problem in its present form. Some

of the useful features of the c" formulation are that its midplane concentration is

bounded even though its surface concentration is zero, that the choice of roots in

equations like eq. (1.31) are easily made, and that it is independent of the species

selected initially to be preserved in rp for the derivation. While these features are not

essential in studying first-order kinetics, they do expedite the analysis for second-order

and third-order reaction systems discussed in Chapters Two and Three.

Effectiveness Factors

Effectiveness factors are commonly used to compare the average rate of reaction

actually achieved to that possible in the absence of concentration gradients. Three types

are encountered: the internal effectiveness factor, r,,, the external effectiveness factor,

fext, and the overall effectiveness factor, overa, = tinotex Only two need be considered

at a time, and these will be Ier,,, and int. The internal effectiveness factor gives the

ratio of the average reaction rate within the slab relative to the reaction rate at the

exposed surface, Lee (1985)

(1L) fr,(x)dx Di, x -
'Iinrt =Lr
rP, L r'Ps


for the slab problem. The dimensionless problem, such as eq. (1.26), is still a slab

problem, with D = L = 1, so


fr"p(z)d 1"z = I
t'int ) r"(z ) 1
r p(z=1) r"p(z=1)

where r"p is the dimensionless rate expression obtained after applying eqs. (1.19) and

(1.20) to rp. Using eq. (1.31) with eq. (1.26) gives

int =

2 2c"(O) c"(0)


for the first-order problem, again showing the need for the negative sign of the square

root to obtain positive, or meaningful, internal effectiveness factors when using c".

Substitution of eq. (1.32) for c"(0) gives the more familiar result

Stanh(O,,) 1.41
(D lst

valid for both reversible and irreversible first-order reactions with the definition of 4,

given in eq. (1.28). For the irreversible zeroth-order reaction discussed earlier the

equivalent result is

TIint = 1 k 1 1.42

Sin- 1 1.43

S= L 1.44

Both the zeroth-order and first-order results share asymptotic behavior of unity for small

Thiele moduli values of ,k or s, and of 1/1k or 1/4w,t for large Thiele moduli values.

This will be discussed in the next section.

The overall effectiveness factor is identical to the internal effectiveness factor for

the case of no external mass transfer resistance. When the Robin's boundary condition

is required, the overall effectiveness factor is defined

rp( c^ ^c...) r (CP=0) 1.45
"overall = is' c int = ext int = ? ", lint

S(dc"/dz) 1.46
overall =
rP(c =c)

The subscript on c"f, the dimensionless free stream concentration of i, is required,

since the relationships used to eliminate j from the reaction rate expression to form r"p

are valid only within the slab and not in the free stream fluid. To use the last part of eq.

(1.45) or to use eq. (1.46) requires more than just considerable care in writing r"p.

The overall effectiveness factor remains species independent, because the numerical

values of the r"p(c"r) are the same for each i, even though the expression for r"p

in the free stream must be modified for each species. This ultimately is due to the

species' concentrations being interrelated through ratios of mass transfer coefficients

instead of ratios of diffusivities between the free stream and the slab surface. This will

be illustrated for the first-order Robin's problem. For the irreversible first-order

problem, note that

Bi,, c" cAdz -tin 1.47


and, since there is no second species to confuse the issue, r"p(c"A) = st c", 1,

both inside and outside the slab, so

Bi.A c f 1.48
Tioverall = 2 v,, 1
stc fA


1 st 1 s 1 1.49
overall BiA BimA C"A BmA 1 int

This sum of resistances format is discussed in Aris, p. 107 (1975). For the reversible

first-order reaction, the expression, r"p(c"fA) = 4t c"f 1, does not express the

dimensionless reaction rate in the free stream nor does r"p(c"m) = st c"C 1,

since both expressions contain no information about the relative rates of external mass

transfer. Both expressions implicitly contain the effect of external mass transfer only for

the selected species. The primary use of eqs. (1.45) and (1.46) in this work will be for

single species reactions. A general expression for reactions with more than one species

requires a new derivation beginning with the definition of the overall effectiveness factor

in terms of a ratio of the dimensional reaction rate expressions, rather than the

dimensionless reaction rate expressions in eq. (1.45). This will be discussed in Chapter


The Generalized Thiele Modulus

In the above section, it was noted that the zeroth-order and first-order internal

effectiveness factors were both asymptotic to unity and 1/4. This was not entirely

accidental. The factor of two in V permits this match. In 1965 the three independent

authors, Aris, Bischoff, and Petersen, developed defining relations for an arbitrary

reaction rate expression that would produce a generalized Thiele modulus, to, such that

,t would be asymptotic to unity and 1/4,,

Lr (e,)
2 f D,(a) rp(a) da

see Lee (1985), where rp has been rendered in terms of a single species through relations

like eq. (1.20) to permit integration. The lower integration limit, .q, is zero for an

irreversible reaction, and a sensibly chosen zero of r,(c) = 0 for a reversible reaction,

usually the one closest to c,. The two moduli derived so far, 4, and 41,, are also

generalized Thiele moduli for the zeroth-order and first-order reaction respectively.

Second-order and third-order power law kinetics for irreversible single species

reactions have been investigated in the context of traditional and generalized Thiele

moduli as part of studies on general n-th order kinetics

Dd2 = ken. 1.51

It should be noted that the stoichiometric coefficient is neglected in writing the right-hand

side. This gives rise to no problem for single species reactions where it can be assumed

to be absorbed into the rate constant, k, although it is a little careless. Because this work

deals with many multiple species reactions, the stoichiometric coefficients will always be

indicated except when discussing earlier work. Making eq. (1.51) dimensionless with

the slab half-width, L, and surface concentration, c,, leads to the traditional moduli

2 L kc= D 1.52

The generalized moduli are very similar in form, needing only an appropriate scaling

term to produce the desired asymptotic behavior

2 (in+1)L2k 1) 1.53

The generalized moduli for n = 0,1,2,3 are given by

2 L2k
Sn = 0 1.54

2 L2k
SG n = 1 1.55

2 3L2kec 1.56
o = n = 2

2 2L2k2, 1.57
= D S n = 3

These four irreversible reaction systems will be used throughout as a basis for

comparison with results obtained later for second-order and third-order reaction systems.

Plots of %, versus to for these values of n date to Bischoff (1965) who determined the

curves for n = 2 and n = 3 by numerical integration of the reaction-diffusion equation.

Figure 3 is a plot of this type generated from analytical solutions derived as part of this

work. The figure is in general agreement with Bischoff's plot as well as one in Aris'

book (1975), but these tend to lose the feature that the curve for n = 2 is closer to the

curve for n = 3 then to the curve for n = 1 over the full range of ,o-.

In the context of the generalized modulus one can say that, at fixed 4 a reaction

giving results similar to those for first-order is less efficient than irreversible zeroth-order

reactions but more efficient than irreversible second-order or third-order power law

reactions. In other words, every physical single reaction system, either irreversible or

reversible, has a location on Figure 3, to and %t, which can be compared to n-th order

behavior. This is perhaps the only setting in which diverse reaction kinetics can be

readily compared. There is one further note concerning Figure 3. It is common to refer

to the region 4,o < 3 as the region of negligible diffusion effects, r, 1, and to the

region 4, > 4 as the region of diffusion control, i.e. insensitive to the rate expression

S e




2 0




* -


form. It is not common, and also not wise, to base these distinctions upon c"(0), since

it may be possible in some systems to have c"(0) small and %, large. The asymptotic

behavior with respect to 4,G is guaranteed, although anything is possible between the

asymptotes. This is amply demonstrated in the literature for n < 0 power law kinetics,

for nonisothermal exothermic reactions, and even certain types of Langmuir-Hinshelwood

kinetics. Results for most Langmuir-Hinshelwood kinetic expressions, however, lie

between the curves for n = 0 and n = 2 in Figure 3.

Shape Factor Normalizations

Zeroth-order and first-order solutions have been derived in cylindrical and

spherical coordinates, as well as some two-dimensional and three-dimensional geometries.

Previous workers have found that a characteristic length choice of Vp/SCx, the particle

volume divided by the external surface area, brings first-order internal effectiveness

factor solutions for slab, cylinder and sphere into approximate agreement, i.e. about as

close as the curves for n = 1,2,3 in Figure 3. There is excellent asymptotic agreement

and approximate agreement between the asymptotes. A more sophisticated

renormalization was proposed by Miller and Lee, (1983). The potential of shape factor

normalization has been somewhat neglected, since analytical solutions existed for zeroth-

order and first-order reactions in many geometries, making shape factor normalization

unnecessary, while there were virtually no other analytical solutions in any geometry.

This work presents a very large number of new analytical solutions in the slab geometry

for which there are no corresponding solutions in any other geometry. An opportunity

exists to extend these new solutions to other geometries at least approximately using the

established principles of shape factor normalization.

Summary and Goals

Ideal catalysts would always have high internal effectiveness factors. A new

catalyst of arbitrary geometry can be compared to an equivalent slab through shape factor

normalization. The internal effectiveness factor of the slab can be calculated

numerically, estimated from a graph such as Figure 3, or calculated from an analytical

solution. Such procedures are facilitated by having a large library of analytical solutions

available, but this has not been the case. Analytical solutions for the internal

effectiveness factor exist for n = 0 and n = 1, as well as one for n = -1 in terms of

Dawson's integral, Aris (1975), plus Thiele's solution for n = 2 power law kinetics

which is

n2D 3L2k o) 1.58

where Co = c(x=0) is obtained from the solution of either

L2 kF sn 3' 41 (s eo)J es 1.59
D 3 F .sin- sin(15)
SO N4 V/-l-e,/0 C

for c,/co < 1 + V3, or

2 I2K(15" ) F sin- o) sin(15*) 1.60
SD 4 1 + /

for c,/c > 1 + V3, where co is the dimensional midplane concentration, F is the elliptic

integral of the first kind, and K is the complete elliptic integral of the first kind. Thiele's

original notation was modified to conform with that used throughout this text. Thiele did

not include an expression for c(x), which presumably would be similar to eqs. (1.59) and

(1.60) with c(x) replacing c0 in one or more places, plus an appearance by x perhaps

replacing L. Eqs. (1.59) and (1.60) are awkward in appearance, require two different

elliptic integrals to evaluate, are each valid over only a specified, but implicit range, and

are algebraically implicit in solution character since c0 is an unknown value of the

solution which is needed in eq. (1.58) to determine t,. The extension of this solution

method to other second-order kinetic systems may be possible, but has not been studied

by others as far as is known by the author. It may well be that this result discouraged

other engineers from attempting to solve more complicated second-order kinetic reaction

and diffusion problems. Nevertheless, one can hardly expect to solve a nonlinear second-

order differential equation without at least having to evaluate one nonlinear algebraic

expression to satisfy the boundary conditions. The following chapters will deal with the

solution of the two problems

d2 = azc2 + a, + a0 1.61


d2c a3 + a2 + al + ao 1.62

derived from net second-order and third-order reversible and irreversible elementary

reaction rate laws using the interrelationships between species developed earlier and using

both the Dirichlet boundary conditions and the Robin's boundary conditions. Analytical

solutions for the concentration profiles will be derived that are valid throughout the

allowable parameter ranges without using elliptic integrals. The solution procedure will

require determination of the midplane concentration in an implicit equation to satisfy one

boundary condition. This single implicit step will complete the concentration profile

equation as well as giving the internal effectiveness factor, which is a simple function of

the Thiele moduli and the midplane concentration. Thus the new solutions will be more

compact and easier to compute than Thiele's solution for n = 2 power law kinetics, as

well as covering a much broader range of kinetic diversity.


Overview of Second-Order Reaction and Diffusion

Limited work exists in the literature on the analysis of isothermal, coupled steady-

state second-order reaction and diffusion problems. Several specific types of kinetics and

physical configurations have been studied, and the results are summarized nicely by Aris

(1975). This chapter will accomplish several interrelated goals. The conventional

second-order reactions, those which are irreversible, reversible with another second-order

reaction, and reversible with a first-order reaction, will be analyzed. The analysis in this

chapter will be limited to the Dirichlet problem for a slab. This will lead to equations

interrelating the individual species' concentrations as well as to an analysis of the

minimum number of parameters needed to formulate a single general problem. This

general reaction diffusion problem will be solved analytically for the dimensionless

concentration profile with the aid of the Weierstrass elliptic Pe-function, p, in terms of

just two parameters. Three different dimensionless concentrations will be used to

facilitate different stages of the development. While this is a little awkward, each has

its advantages and disadvantages. These will be discussed. The mathematical nature of

the Weierstrass elliptic Pe-function will be outlined. The general nonlinear governing

differential equation leads to an implicit algebraic expression to be solved to satisfy one


boundary condition. The dependence of the solution of this equation on the two

parameters will be examined. Second-order reactions reversible with another second-

order reaction possess a special feature. Under certain values of the kinetic rate

constants and diffusivities, the governing differential equation will degenerate to first-

order kinetics with linear differential equations. This could be forced to occur in many

instances by adjusting temperature, for example, without invoking any infinite

diffusivities or other asymptotic limits. Finally, an analytical expression for the internal

effectiveness factor will be examined as a function of the two most fundamental


Burghardt's recent article (1990) still claims that there is a lack of an explicit

relationship between the internal effectiveness factor of a pellet and the characteristic

parameters for other than zeroth-order and first-order reactions. Burghardt uses theorems

from singularity theory to divide the feasible parameter space into regions with different

bifurcation diagrams, i.e. with different numbers of possible solutions, for an exothermic

reversible first-order reaction with a uniform pellet temperature different from that in the

surrounding fluid. Their analysis should now be possible for selected second-order

kinetic cases. Relationships between the internal effectiveness factor and the new

analytical solution will be derived here for seven classes of second-order kinetics.

Analytical solutions for second-order kinetic expressions are generally only possible for

Cartesian coordinates and then only in terms of elliptic functions or integrals. These are

either sufficiently unusual or unavailable that many engineers would undertake a

numerical solution of the problem instead. This need not be the case as will be shown


below. Furthermore, a calculated concentration profile is not required to determine the

internal effectiveness factor. The internal effectiveness factor can be calculated if the

midplane concentration alone is known.

As will be seen below, the initial dimensionless forms of the individual species'

governing equations need contain no more parameters than the number of distinct species

in the elementary reaction plus one for the Dirichlet problem, i.e. up to five for the most

general second-order reversible case, A + B 4 P + R. The number of dimensionless

groups contracts to just three species' dependent parameters in the first step of the new

analytical solution procedure for the dimensionless concentration profile (fewer for

irreversible reaction cases).

The analytical solutions derived in this chapter contain an expansion in terms of

the Weierstrass elliptic Pe-function that is arguably easier to use than those for the

traditionally acceptable first-order irreversible reaction in an infinite cylinder, which

requires the modified Bessel function of the first kind.

Some earlier analyses exist in the literature for simultaneous second-order

reactions and diffusion, Thiele (1939), Maymo (1966), Bailey (1971), Tartarelli (1970).

The previous analyses have not produced much general theory about second-order

elementary reactions. These early analyses have tended to be in the spherical coordinate

system of catalyst pellets and have relied heavily on numerical methods to obtain

solutions. Often these authors have invoked stoichiometric surface concentrations and/or

a common diffusivity for all species to reduce the dimensionality of the parameter space

in order to permit graphical presentation of the results. Indeed, one recurring problem


in previous presentations of solution results for simultaneous diffusion and second-order

kinetics has been the confusion in identifying the minimum dimension of its parameter

space. Among these earlier analyses are some with up to eight dimensional parameter

spaces. The first step reduction of the governing equations taken in this chapter produces

a three-parameter basis for isothermal reaction kinetics with second-order forward

reactions, and either second-, first-, zeroth-order or no reverse reactions (nine kinetic

schemes are possible, counting two with zeroth-order reverse reactions). Numerical

values of the three parameters depend on the species selected for the calculation. This

is not a major drawback, but it can be inconvenient. All three parameters are nonzero

except when, (1) a special first-order degeneracy occurs, (2) the surface concentrations

reflect complete or equilibrium conversion, or (3) the reaction is irreversible (2A -*

products and A + B products).

Further analysis of this problem yields a two-parameter basis which is independent

of the chemical species selected for the derivation. These two parameters are considered

to contain the essence of second-order reaction and diffusion behavior. A clear picture

has emerged with respect to the conditions under which the general second-order

reversible reaction system is free from diffusion limitations or is diffusion limited as a

result of the reduction to two parameters. Also emerging later in this chapter will be

some less obvious results. These include the manner in which the internal effectiveness

factor is affected by the two parameters, bounds on the dimensionless midplane

concentration and its estimation, a bound on the relative magnitudes of the two


parameters and asymptotic behavior of the internal effectiveness factor as a function of

the two parameters.

Some additional practical problems have been studied by others. These will be

briefly summarized here. Bailey (1971) has analyzed the problem of optimum surface

concentration for the irreversible reaction case A + B -- products, when the diffusivities

are unequal. Tarterelli et al. (1970) have studied the two irreversible kinetic cases in a

macropore-micropore catalyst model. Thiele (1939) developed an analytic solution for

the slab midplane concentration and internal effectiveness factor for the single species'

irreversible reaction, 2A -- products, using elliptic integrals instead of elliptic functions;

see Chapter One. The selectivity of porous catalysts undergoing parallel reactions was

examined by Roberts (1972), and includes second-order reactions in parallel with either

a zeroth-order or first-order reaction. Roberts also used elliptic integrals. Two early

studies on nonisothermal behavior are noteworthy in that the authors disagree on the

results obtained, Tinkler and Metzler (1961) and Petersen (1962). These earlier studies

as well as this study neglect the effect of volume change during reaction, which is most

often valid when the system concentrations remain near those at the surface, when the

effective diffusivity is predominantly due to Knudsen diffusion, or when the fluids are

dilute gases in an inert gas solvent or liquids.

Governing Equations and Dimensionless Groups

Consider the general second-order reversible isothermal reaction of the type A +

B # P + R, occurring in a slab-like medium or along the axis of a cylinder as


illustrated in Chapter One, under the assumption of constant effective diffusivity, Di.

The governing equations take the following form

d2 A d 2p
DA --d -r, DP r,
dx2 dx2

DB rp DR rp
dx2 dx2

where rp = k,cpCR lkACB and x is measured out from the slab midplane. These

expressions can be written more generally as

Si r 2.2

The rate of product formation rate expressions, rp, considered in this chapter are given

in Table 2.1. The first seven stoichiometries listed there are discussed in this chapter.

Zeroth-order reverse reactions are not discussed here. The analysis that follows always

assumes a symmetric problem, i.e. the slab with sealed midplane or the slab with two

exposed faces with identical surface conditions. This leads to either a Dirichlet or

Robin's problem at the exposed surface. The Dirichlet problem is the more basic of the

two. Generally a solution can be found for the corresponding Robin's problem, if the

Dirichlet problem can be solved.

There are many steps between this point and the analytical solution. The solution

strategy to be employed will be outlined here so the reader knows what to expect in the

pages ahead. The initial step will be to take a set of n coupled nonlinear ordinary

Table 2.1 Typical Second-Order Kinetic Schemes












A+B P + R

A+A P + R




A + B-

A + A-

A + B zeroth-order

A + A < zeroth-order

Rate Expression

rp = k1CCR IACB

rp = kcR kfcl

rp = kCp- ktCACB

rp=kjCp-k ktC
rp = k ACB

rp =- k K

rp = k, kfcAB

rp = k, kt
r = -


differential equations (ODE's) and form n single species ODE's, where n is the number

of species appearing in the rate expression. These uncoupled ODE's will be made

dimensionless using a characteristic length and the species' surface concentrations. It

will be prudent at this point to condense and label numerous common dimensionless

groups contained in the ODE's to simplify notation. The fourth step is an optional step

as far as obtaining the analytical solution is concerned. A properly chosen linear

transformation of the dimensionless species concentrations leads to a species independent

ODE containing fewer parameters which are themselves independent of species in the

sense that they are identical regardless of which species is used as a starting point for the

transformation. The variables are separated and integrated one time by conventional

methods in the fifth step, see Lee (1985), leading to a concentration polynomial one

degree higher than in the original reaction rate expression. Sixth, take a complete, or

entire, Taylor series expansion of this polynomial about one of its zeros (no truncations).

Next, make a very special and particular nonlinear transformation of the concentration

variable. Then separate variables for the second time and recognize that the result is

something with particular significance, i.e. the Weierstrass elliptic Pe-function. Finally,

back substitute to whatever point desired to produce dimensionless or dimensional

concentration profile equations.

This chapter will consider the Dirichlet problem. Conventional boundary

conditions for the Dirichlet problem are

ci = ci, x = L

d = 0 x = 0. 2.4

The following relationships between the species are readily derivable from the governing

equations, see Chapter One, since (Dici/vi Djcj/vj) are solutions of a potential function,

regardless of the coordinate system chosen, and thus equal to a constant. For the

Dirichlet problem, this constant is known from the surface condition, i.e.


DAA + D, = DA + D, e 2.6

and similarly for R as for P. These can also be written more generally as

D. D.
--(c- c) 2.7
v. vj

These expressions are critical to continuing the analysis for nonlinear kinetic rate

equations. All but one concentration can be eliminated from rp using eq. (2.7). Then

the system can be reduced to a single ordinary differential equation. If A is selected, and

concentration and length are made dimensionless with CAs and L, then the above second-

order reversible system reduces to

d2cA (L2kA L2krDAC As 2 L L2k

kes + + R ,L A Ak R a 2.8



where CA = C^/C, and z = x/L. In spite of the lengthy coefficients, eq. (2.8) is an

improvement over eq. (2.1). The right-hand side of this equation will be referred to as

the diffusion modified reaction rate expression. Dimensionless species' solutions depend

only on the values of the three coefficients of cA, CA and 1. The above equation can be


d 6a c,2 + 6acA + 2a 2.9

with revised boundary conditions

CA =1 z= 1 2.10

0 z = 0. 2.11

The numerical constants in eq. (2.9) are per Whittaker and Watson (1950) and, when

following their solution method, these choices reduce the occurence of fractions in the

Table 2.2 Coefficients of c2 for Second-Order Reactions

Case Species 6ai

S+^ 2 t2 t2
SB B+- rApB +A +R/ A+B

I P #+- A+B/rAP

I R ~, >2 2PAP+

II A 2(2A rAPp+R/4)

II P I+R- 42A/rA

II R R+I 4A +PrA +R

III A 2(i2A AAP 2P)

III P 2(2&p 42A/rAP)



IV P -$A+B/rAP

V A 22A

V P -442A/rAP



VII A 22

Table 2.3 Coefficients of c, for Second-Order Reactions

Case Species 6

I A I+A- +B + P+R + R+P+2FAPp+R

I B $+B-B+A +"P+R +R+P+2AP B+AR +P/A+B
ID 2 t +2 2 +
I P R+P-4 +R +A+B +B+A+24+B/rAP

I R Pi+R- R+P+ A+B B+A+ B+2 AR+P/AP P+R

II A P 2+R + p+rAP +P

II P D2 2
II P RLP-IP+R+4 2A 82A/rAP

II R IP+R-,pR + 42A + 82 +P/ FApP+R

III A 2(242p + 2rFApp)

III P 2(24 +24A/AP)

IV A 2 +A B+

IV B A+B- B+A+
IV P (D 2 +D2 t2
IV P 4P+A+B' B+A+24 A+B/PAP

V A p

V P 4+4 +4 A/rAP


VI B 2A+B- +A


Table 2.4 Constant Coefficients for Second-Order Reactions

Case Species 2a,

I A -(4+p++ R+rAP P+R+ R+p/'AP)

I B -(42+R+ 2+PrA +A 2 +R/4 +B+ +B +P/rAP2+A

I P -(IA+B B+ArAP B+A ++B/rAP)

I R -( B+ B)+A+ rAPt+A P+R/,R+P+A +B R+P/rAPDP+R)

II A -(p+R + +P+2p +p/FAP+ Ap4p+R/2)

II P -_A(4+PAp+4/PAP)

II R -4A2(4 +rAPp +R/42 ++4 +p/rPAP +R)

III A -2 4p(2+rAP+ /rAp)

III P -2i2A(2+PAP+ I/FAP)

IV A (1 + 1/PAP)

IV B -4p(1 +4AB/FAP B+A)
IV P 2 4p2 pot2 _4p2
IV P 4A+ B +AAB+A A +B/rAP

V A -4p(1+/rAp)

V P -2A(4+rAP+4/rAP)

VI A 0

VI B 0



final results. Six dimensionless groups appear naturally in these three coefficients in eq.

(2.8). If B, P or R had been selected instead of A, then some of these groups and some

new, but similar, groups would have appeared. Expressions for the ag are given by

reaction and species in Tables 2.2-2.4 for the seven cases free of zeroth-order reactions

using a notational standard described below.

A shorthand notation system would be useful to help minimize the number of

dimensionless subgroups appearing within the more general parameters and facilitate

studies of parameter interrelationships. After all, there are 17 species dependent ODE's

under consideration at this point. Such a notational system should also retain some

physical significance or familiarity. Toward this end, a set of irreversible reaction Thiele

moduli will be defined

I2 = L2 ki,/ D 2.12
I i LJ

where i,j are either both products or both reactants. The kl indicate the rate constant,

either kf or k,, multiplying the concentration of species I in the rate expression. Such a

group would appear naturally in the study of the irreversible reaction I + J products.

Large values of such a group would suggest possible strong diffusional effects, while

small values suggest a possibility of nearly diffusion free operation. The true state of

affairs will be revealed by further analysis. If i = j, define the irreversible reaction

Thiele modulus

0 = L2 ki ,/D, 2.13

which would arise in studying 21 products or I + I products, and define

r Di is 2.14

where i and j are not both products or both reactants. This ratio serves as a bridge

between product and reactant species. The set {4A+,, IB+A, P+R, IR+P, IAp} is

sufficient, as can be seen from the tables, to write the dimensionless governing equations

of all four species in the reaction A + B ; P + R. Another sufficient set {4A+B, 4,+R,

"AB, rAP, IAR} is prototypical of many similar sets, as is {lk/k,, 2A+B, AB, AP, ,AR }

The first set is also sufficient to express the internal effectiveness factor using the

derivative of the dimensionless concentration. It can be seen in Table 2.2 that 4'+j is the

coefficient of ci when the rate constant of the opposing reaction is zero, and that it has

similar dependence on physical parameters as the generalized Thiele modulus for a

second-order irreversible reaction of a single species given in Bischoff (1965), Aris

(1975), and Lee (1985)

2 = 3L2ke,/2D. 2.15

The following interrelationships between the parameters in the new notational system

should be noted

Fi = 1/r.. 2.16

r = -j k 2.17
= i


For the case of mixed first-order and second-order reversible kinetics, a first-

order irreversible reaction Thiele modulus will be defined in the traditional way as in

Chapter One, ,I = L2'k/Di, where i indicates the species undergoing a first-order

reaction and its rate constant. With this adopted nomenclature, the coefficients of the

governing equations for many second-order systems can be readily condensed.

Expansions for the coefficients, aj, for the first seven second-order kinetic expressions

listed in Table 2.1, have been formulated in terms of the irreversible reaction Thiele

moduli and a single FIj. The results of employing the above nomenclature scheme are

summarized in Tables 2.2-2.4, for d2ci/dz2 = 6a1,ic + 6a2ic, + 2a3,.

Second-order systems with a zeroth-order reverse reaction can be handled by the

integration methods that follow, but are not considered in this chapter. In dealing with

such reactions, it is potentially possible to associate a real species with the zeroth-order

reaction, and to have the concentration of that species drop to zero at an intermediate

point in the slab. This requires the introduction of an additional parameter to establish

this position, and conventionally the rate constant for the zeroth-order step would have

to be made zero or set equal to the opposing rate beyond that point of penetration.

Solution of the Dirichlet Problem

The 17 dimensionless species dependent governing equations in the above seven

conventional cases are of the same form, namely

dc -= 6 aic + 6a2,c + 2a3, 2.18

where ali, a2i and a3, are those coefficients summarized in Tables 2.2-2.4.

Two generic ways are developed here to further unify the species dependent

parameter problem of eq. (2.18) to a unified two-parameter boundary value expression

for the general case where the aij are all nonzero. This step is not critical to solving the

differential equation, but reducing the parameter space dimension is quite useful and well

worth pursuing. For the first of these methods, let

c' = 6a,,(c, 1). 2.19

The variable, c', could also be interpreted as a dimensionless extent of diffusion reaction

variable, since it characterizes the species independent evolution of the diffusion reaction

and has the value of zero at the slab surface. The result that c' is independent of species

is anticipated in leaving it unsubscripted. Then eq. (2.18) becomes

d2 = c2 + (12ali+6aZi)c' + 6ali(6ai,+6a2+2a3,) 2.20

which can be rewritten as

d2c = 2 + 0 + T 2.21

with boundary conditions

c' = 0 at z = 1 2.22

d 0 at z = 0 2.23


where ,t = 12a1, + 6a2i and = 6a1i(6a1i + 6a2i + 2a3). The sign of c' is no longer

necessarily positive, unlike the case for c, in eq. (2.18). The scaling choice made in eq.

(2.19) leads to a coefficient of unity for c'2. Other scaling choices are possible that

would lead to a fixed coeffient other than unity for c'2, but they do not seem to offer any

particular advantages.

The procedure for the first integration of the general ordinary differential equation

above is well known; see Lee (1985). Let q = dc'/dz, then

dq c'2 + s2c + 2.24
d z .

Multiply both sides by 2qdz, equivalent to 2dc', and integrate to obtain

q2 = 2 '3 + ,c + 2' c' + a 2.25

where a is an undetermined constant of integration. It is possible to substitute the

definition of q, separate variables, choose integration limits, etc. at this point, but the

analytical solution of that integral would not be known. Note that q(c'(0)) = 0 by the

boundary conditions. This establishes a relationship between a and c'(0), namely

2 22.26
a = --c'(O)3 ~tc'()2 2Yc'(0) 2.26

The value of q(l) will be needed later for determining the internal effectiveness factor

which depends on the surface concentration gradient. It is given by

q(1) = i = -c'(0) 0 c '(0)2 2c'(0) 2.27


where the same sign as is generally correct physically. If a reliable estimate of c'(0)

is available, and concentration profiles are not required, then it is possible to avoid

performing the second integration, and to proceed directly to the calculation of the

internal effectiveness factor. Assuming the contrary, it is important to recall that c'(0)

is a root of q2(c'(0)), since (c'" c'(0)) = (c' c'(0)) times a (polynomial of order

n-1) for all relevant n in this study, CRC Standard Mathematical Tables, (1975). The

fact that the remaining boundary condition is at z = 1, not at z = 0, must be temporarily

ignored in the following steps. First expand q2(c') about q2(c'(0)) in a Taylor series.

q2(c') = (2c'(0)2 +2~ 0f c'(0) + 2')(c' -c'(0)) +
(4c'(0) + 2 ,)(c' c'(0))2 + 4(c' -c'(0))

Note that this equation is complete, and not an approximation, since eq. (2.25) is a cubic

polynomial. The point of this step is to eliminate, or at least disguise, the constant term

with c'(0) in eq. (2.25) by moving it into the differences with c'. Then define

q2(c') = 4A3(c' c'(0)) + 6A2(c' c'(0))2 + 4A1(c' c'(0))3 2.29

Substitute dc'/dz for q, separate the variables c' and z, and then take the integral from

(0, c'(0)) to (z, c').

S= fda 2.30
0 c0) )4A3(a-c'(0)) + 6A2(a-c'(0))2 + 4A(a-c'(0))3

Make a change of variables to r = (a c'(0))':

z= dr 2.31
t 14A33 + 6A22 + 4A,

Note that if the original kinetics had been third-order in c', it would only have added

a constant term to the quantity in brackets, which would have been called Ao; see

Chapter Three. This is therefore a powerful substitution, which is not fully exploited

here but will be in Chapter Three. The next substitution is conventional and

simultaneously reduces the coefficient of the cubic term to four and eliminates the

quadratic term. Let a = A3r + A2/2, then

z = do 2.32
o 403 (3A2 4A1A) a (2AIA2A -A3)

Sf do. 2.33
a 44o3 -g2 -g3

The numbers, g2 and g3, are invariants of the original third-order polynomial, Whittaker

and Watson (1950), and are also (and more readily) given by

( s,-4 Y) 2.34

4c'(0)3 + 6 t c'(0)2 + 12 c'(0) + 6 ,s't 6t 2.35
3 216


The relationship between a and c'(0) is needed to verify the g3 equivalence. The

invariant, g2, is independent of c'(0) and therefore ab initio known.

Eq. (2.33) is one defining relation of the Weierstrass elliptic Pe-function. It is

also written with an elaborate script capital P, p, Whittaker and Watson (1950). The

variable transformations following eq. (2.30) were made with this end in mind. The

cubic polynomial in the denominator will be referred to as the canonical polynomial of

the problem, even though it is not quite the canonical form of a cubic polynomial given

in most mathematical handbooks. These lack the factor of four and use addition rather

than subtraction of the lower order terms. The defining relation for the Weierstrass

elliptic Pe-function translates eq. (2.33) to

S= p(z; g2,g3) 2.36

Substituting r for a and then c' for 7 gives the analytical solution for the dimensionless

concentration profile in terms of the still undetermined midplane concentration

6(c'(0)2 + l2,,c(0) + )
c'(z) = c'(0) + 6( 1c'() + 2.37
12 P(z; g2,g3) 2c'(0) 4I,,

where c'(0) must be evaluated from eq. (2.37) using the boundary condition of eq.

(2.22). Eq. (2.37) can be inserted into eq. (2.19) to give the dimensionless species

concentration profiles.

An alternative approach, Chapter One, to reducing eq. (2.18) is to take

c" = lc c' 2.38
6al, + 6a2, + 2a3, T

or equivalently in terms of original problem parameters


C, ( i ^)
L2' vi r

where rp, is rp evaluated at the surface concentrations. This leads to the following

revised governing equation

d2c" cu+c"
dz2 c' + 1St c 1


where 1,, and I are identical with those defined by eqs. (2.20) and (2.21). The

boundary conditions for c" are identical with those for c'. The analytical solution for

eq. (2.40) is derived exactly as eq. (2.37) and is given by

c"(z) = c "(0) +

6(- c"(0) + c "(0)- 1)
12p(z;g2,g3) + 2Yc"(0) 2,

with invariants

4lst 4 T)
82- -- ^ -
9g2 12

-43c"(0)3 +6V2 tc" (0)2 12 2c"(0) + 6T2t lst
3 = 216





The same results could be obtained by direct substitution of eq. (2.38) into eqs. (2.34),

(2.35) and (2.37). The numerical values of the two invariants are unaffected by the

change of variables. The integration technique could have been applied to the original

ODE with three parameters as well. Remember that the reduction to two parameters was

not required to solve the ODE. Interested readers may consult Koopman and Lee, 1991,

where this form is derived in full. Inclusion of this integration is redundant in this

chapter, and the three-parameter form is also less useful for what follows with regard to

second-order reaction and diffusion systems than the two-parameter form. Using either

eq. (2.37) or eq. (2.41), it is now possible to present general solutions for internal

effectiveness factors, slab midplane concentrations, etc. as functions of 41t and \I.

The Two Fundamental Parameters

The two fundamental parameters, ,t and *', have significance beyond reducing

the problem from one with three parameters, the aj, to one with two parameters. The

first of these parameters, 4{, is never negative and is numerically independent of the

species chosen in the derivation. As seen in Table 2.5, it takes the form of a sum of

irreversible reaction Thiele moduli, i.e. moduli that one would derive from each species'

governing equation by ignoring reverse reactions

s, = (vi ri)2 2.44

where the sum is over all species, and ,o i is either 4,+j, '2i, or 4,i depending on the

kinetics. This expression is reminiscent of the one for the Thiele modulus of a first-order


Table 2.5 Defining Relations for the Parameters of the Governing Equations for Seven
Second-Order Reaction Systems


,2 it2 +iA2 +4h2
A+B+ -+B+A P+R -R+P

42A +P+R + +P

4.tA + 404p

2 2
SA+B + IB+A +P

4 2A + 42










kA+BB +A +P+R4R+A-rAP~B+A P+R- A+BR +P AP

42A +p RR+p-rA 2A P+R-4A I'+P/rAP
44t4+4.t4 -42 2

4 +442, 4A2P(rIA+ I/PAP)

2 t2

4A 4A2A AP


4 .A


reversible reaction, 4c = -A + Cl given in Satterfield (1970) and in eq. (1.28). A

Thiele modulus symbol has been chosen for ist because of the direct relationship with

irreversible reaction Thiele moduli in eq. (2.44), as well as because it multiplies the first-

power concentration term. It also has a similar influence on calculated solutions as that

of the Thiele modulus in first-order problems when higher order effects are minor. The

first-order modulus, 4,, is also independent of the scaling choice of eqs. (2.19) or

(2.38). A general derivation for the four species reaction, vAA + vBB # vpP + vRR,

and rp = kl C~ R kCACB, gives the following expression for V,,

2 VpL2kdp, + vRLk2R, e vAL2kA vBL2kfB 2.45

Simple rules apply for deriving the expressions in Table 2.5 from eq. (2.45) when, for

example, A and B are identical species, or when the reverse reaction is first-order (no

R). Since complete expansions are given in Table 2.5, discussions of these rules are

deferred until Chapter Three and the more complicated expressions for third-order


The second dimensionless parameter, 6a1,(6ai, + 6a2i + 2a3) or NI, has no first-

order analog, and in eq. (2.40) it is the coefficient of c"2. Symbolic expansions for *I

are given in Table 2.5. There it is seen to be composed of irreversible reaction Thiele

moduli as t,, plus the bridging parameter, lAp. The irreversible reaction Thiele moduli

appear only as products with themselves or other moduli. For the general second-order


problem just discussed above, can be expressed in terms of original physical constants


P = L VPVRkr VAk )(CR -kfces)e 2.46
Dp R DAD, )rAe 2.46

It is seen that is composed of a difference times the reaction rate at the slab surface.

Both parts of this expression are of arbitrary sign, so can be either positive or

negative. In fact could be zero for nonzero surface reaction rates through the first

term in eq. (2.46). This leads to precise first-order behavior of the solution as will be

discussed later in this chapter. Note that for irreversible reactions, k, or kf is zero, and

* is always positive. The parameter, *, will be referred to as the second-order

modulus. It is dimensionally similar to the fourth power of a Thiele modulus, L4k2',/D2.

The second-order modulus, *, is independent of the chemical species selected as a

derivational basis as was the case for V,, and unlike the case for the numerous aij. A

discussion of simple rules to use eq. (2.46) to form the expressions in Table 2.5 will be

deferred to Chapter Three as mentioned above.

The feasible ranges for 4,, and define a region somewhat smaller than a semi-

infinite plane, (4it 2 0, f). Although can be positive or negative, its positive range

is limited physically by the condition that all species' concentrations, ci(x), must be real

if the diffusion modified reaction expression (right-hand side of the diffusion equation)

is set equal to zero. The dimensional, diffusion modified "equilibrium" concentrations

are related to the dimensionless diffusion modified equilibrium concentration, c",(z),

which constrains c"(z) by

2 Y
-st + Ols -. 4 2.47
0 c c"(z) < c"eq(z) = 2-2.47

where the right-hand side arises from applying the quadratic formula to the right-hand

side of eq. (2.40). This diffusion and reaction driven concentration shift generally does

not lead to an equilibrium state if applied to the set of external surface concentrations in

the absence of diffusional effects. In other words the extent of reaction of the diffusion

modified rate expression is a physically different concept from that in homogeneous

kinetics. Applying the requirement of real valued concentrations at equilibrium leads to

the following physical parameter constraint

Slst 2.48

which ensures c",(z) is real and also implies g2 > 0. The same result, i.e. that ,I/4

is the largest value available to I, is also obtained by solving the problem: maximize

*(4,+B, I,+A, Pp+R, R+P, rAP), subject to constant ,,, and the species' parameters

constrained to be positive.

An additional word is in order about the constraint in eq. (2.48). The equality

in eq. (2.48) can be examined by expanding {,,/4 I for Case I of Table 2.1 in the

physical constants of the original problem to give

S-4 = L'4 k ( e _( eft) + k,( Ps ) +
1st D D DR D
D ADB'2.49
L'kfk (DA A + DP P) (DcBS" + DR ) > 0

If species A and B are identical and k, = 0, or if cA,/D = CB,/DA and k, = 0, then 14/4

- 9 = 0. Similar conditions occur if the roles of (A, B, k,) are switched with (P, R, kf).

Therefore the equality of eq. (2.48) always holds for the case of the irreversible reaction,

2A -, products, and can occur in the case of the irreversible reaction A + B -- products.

The equality is not expected to ever occur for reversible reactions. This reflects the fact

that even if products are not present at the surface, they will be present in the slab and

will influence the achievable reaction rates in the slab.

The slab midplane concentration, either c'(0) or c"(0), must be obtained by

solving a nonlinear equation to complete the analytical solution of eqs. (2.37) or (2.41).

The development so far is such that c"(0) has a much narrower range than c'(0), which

ranges from +oo to -oo in the ( *,t, I) parameter space. The sign of c'(0) changes

with the sign of I, and when I I is large, c'(0) | is also large. The rescaling to c"

in eq. (2.38) yielded some surprising results as seen in Figure 2.1, where only feasible

combinations of 4,,, and I were used per eq. (2.48). Not only was c"(0) always

positive, as expected, but it was seemingly bounded above by one-half. Integration of

eq. (2.40) assuming the right-hand side is a constant equal to -1 gives one-half as the

midplane concentration. If the magnitude of the dimensionless reaction rate decreases





~J 1)
0 0 c;0


9 I-

/ '~a-



U? ~ ~ 6--

o -





s I





o* I
6i u

Ln ~*S







1 0

no 4


in the slab, the midplane concentration falls below one-half. Figures 2.2 and 2.3 provide

some two-dimensional perspective on the surface of Figure 2.1 for selected parameter

values. The limit of 41,t 0 at constant 9 gives values of c"(0) that are within 0.5%

of those shown for 4, = 0.01. This limit has more relevance in less conventional

kinetic schemes.

Figure 2.1 is somewhat deceptive. Since the surface concentration, c"(1),

equals zero, the maximum difference between c"(l) and c"(0) occurs at parameter

ranges that correspond to high internal effectiveness factors. This is counter-intuitive in

the sense that small concentration differences from surface to center are normally

associated with diffusion free cases and large concentration differences are normally

associated with diffusion limited cases. This phenomenon has to do with the inverse I

scaling constant in the definition of c".

Eq. (2.47) had given one upper bound for c"(0), c"q(z=0), but it is only

closely approached when c"(0) goes to zero. For example, at (c1,, ') = (0.1, -105),

c",(0) = 91.6, implying 0 < c"(0) < 91.6, which is a much broader range than

zero to one half. Figure 2.4 plots the scaled departure of either c'(0) or c"(0) from

the diffusion modified equilibrium conversion as given in eq. (2.47) and where

,C-'st -t -4 2.50
eq 2

for c',. Note that (c'e c'(0))/c'q = (c", c"(0))/c" The departure from

equilibrium conversion is negligible for strongly diffusion-limited cases. It is not


formally proper to discuss the solution of eq. (2.40) at I = 0, since the equation was

derived using division by *. The solution to eq. (2.21) for c', however, passes through

c'(0) = 0 at = 0, changing sign in the process, i.e. c" approaches zero over zero

as I -- 0.

Figure 2.4 further substantiates one of the claims made earlier on the behavior of

the slab midplane concentration. The surface shown characterizes the approach to either

diffusion modified equilibrium conversion (for reversible reactions), or complete

conversion (for irreversible reactions), at the slab midplane. This approach is quite

strong as 9 -oo and/or 4k,, -* oo. When -, t,/4 (g2 = 0) as 4,, oo; however,

the approach to high midplane conversions appears to be retarded, i.e. along the diagonal

boundary of the parameter space. The accuracy of eq. (2.54), below, for ci(0) is also

qualitatively indicated in Figure 2.4.

If one defines

1 1 / cosh(Os,)
c" = lim(c"(0)) = 2.51
-0 Ol

then the solution of eq. (2.41) for c"(O) should match that of the first-order problem,

see Chapter One, for

d 2 ,stc 1 2.52

given by c"' above. The group, (co' c"(0))/cc', represents a scaled measure of the

departure of the dimensionless midplane concentration for second-order reactions from


C.o c












that of the first-order problem. This group is plotted in Figure 2.5 against 1,, and *.

The plot includes contours showing the demarcation between values within 1% of first-

order behavior and beyond. Figures 2.4 and 2.5 are somewhat complimentary in their

regions of small departure. Eq. (2.52) gives c"(0) as a function of i,, at = + 10-

to within 0.1%. The double limit II| -* 0, I,|., -* 0 gives c"(0) = 0.50 exactly.

Numerically, c"(0) < 0.5 at all 1092 values of 4,, and examined.

Combining eqs. (2.19) and (2.51) gives an approximate relationship between the

dimensionless species' concentration, c,(0) and ai,, 1,,, and *

F( 1-1/cosh(i ))
4,(0) c,, 1- .--- ) 2.53
6i 6li, ist

approximately valid for |I < 0.1 by derivation, but applicable at some larger \II\

when 4),t is also large per the limits in Figure 2.5. Alternatively, when (c"q -

c"(0))/c"q in Figure 2.4 is small, c"(0) is well approximated by eq. (2.51). This

leads to a second approximate expression for ci(0)

12 ao.
( 2 /
ei(0) 1 eis I + 0 12 aJ 2.54

which is acceptable for it,, > 5 or < -1000. Eqs. (2.53) and (2.54) together cover

much of the parameter space. Therefore, if one does not wish to work with the

analytical solution, then eqs. (2.53) and (2.54) offer one alternative method for

estimating the midplane concentration and ultimately the internal effectiveness factor.












e, 0



0 (0 CM




Finally, one should not generally think in terms of ,,st going to zero at constant

'. This can be seen by noting in Table 2.5 that when < 0, | <- M s,,t, where

M is a constant, such as rAP/2 for Case I (when rAp > 1). In other words, the

magnitude of is generally constrained by the magnitude of t'I and can not be chosen

independently in most actual physical systems; however, there are possibilities for

exceptions. The point is fairly subtle. How can 9 go to zero much more slowly than

?',? First note that 4,s goes to zero if (1) L 0, (2) kf -* 0 and k, 0, (3) Di oo for

all species, or (4) CA,/DB 0, CB,/DA 0, p,/DR 0, and CRs/Dp -" 0. Then note that

I goes to zero under the first three limits above, but the fourth must be modified to (4')

CAsCBs/DADB 0, CpRs/DDR -0 0, CAsBs/DpDR 0, and CPCRs/DADB 0. Since only the

first half of condition (4') is forced by condition (4) for 1,t, it is still possible for I I

to remain large, but the rest of condition (4') often follows from condition (4) in

practice. Still the key terms are CA,CB,/DpDR and CP,CRs/DADB. Consequently, only under

unusual circumstances will I be large when 41,, is small.

The Weierstrass Pe-Function. oi

The Weierstrass elliptic Pe-function, p, is discussed in Whittaker and Watson

(1950) and other texts. One of its original uses was in the study of dynamics of

nonlinear pendulum motion, with friction dependent on displacement squared and/or

cubed, and the sin(x) expansion increased to include the cubic as well as linear term.

In this situation the Weierstrass elliptic function provided an analytical solution to an

approximate physical model. Here, if the assumptions on diffusivity, etc. are valid, the


Weierstrass elliptic function is an analytical solution to a true physical system, since the

reaction rate forcing function is not approximated. The Pe-function, p, is a doubly

periodic function in the complex plane, i.e. p(z + 2mc + 2nwo) = p(z) for any

integers m,n. The ratio of w /ow2, the two half-periods, is not purely real, and the only

singularities in the finite part of the plane are double poles at z 2mw, 2nw2 = 0. The

two periods, 2 0 and 2c2, define a single function. Many equivalent period pairs can

define the same function; however, all such pairs are linked by common invariants, g2

and g3. Abramowitz and Stegun (1964) give a procedure for determining a base pair of

frequencies from the two invariants, but this is not germane to the present problem.

They also give an expression for p in terms of the two invariants, g2 and g3,

1 2t-2 2.55
P(z; g2, g3) = 2 + 2.55
Z k=2

where there is no constant term, and

C2_ 2 C g 2.56
20 28

3 E Cm C-m 2.57
Ck= m=2 k 4
S(2k + 1)(k 3)

Eq. (2.55) is the definition of choice for evaluating the Pe-function in this work.

Abramowitz and Stegun (1964) give the coefficients up to C,9 in terms of C2 and C3 only.

If g2 = 20 and g3 = 28, then C,9 = 6.04*10-6. (Evaluation of any midplane

Values of the Weierstrass Pe-Function


g2 g3

0.0 1.0

0.0 3.0

1.0 1.0

1.0 -1.0

-1.0 1.0

-1.0 -1.0

0.1 0.1

-0.1 -0.1

3.0 -3.0

5.0 2.0

5.0 -2.0

20.0 28.0

20.0 -28.0

40.0 -20.0
















First Five

Terms of p































Table 2.6 Typical


concentration from the final boundary condition requires summing the series at z = 1.)

Some representative values of p(l;g2,g3) are given in Table 2.6. This summation will

be time consuming, except for a digital computer, if g, and/or g, are large in magnitude.

Note, however, that if (g2g3) < 0 and Ig2/g3 | = 0(1), then p(1;g2,g3) and p(1;-g2,-g3)

are closer to one than p(1;g2,-g3) and p(1;-g2,g3). In many regions of engineering

interest g2 > 0 and g3 < 0. Convergence of the series for 0 < z 0.8 can be quite

rapid, regardless of the invariant magnitudes, because powers of z tend to dominate over

products of the invariants.

The two invariants, g2 and g3, are approximately of the same order of magnitude

as the largest irreversible reaction Thiele modulus, e.g. 4+A,, raised to the power of two

or three respectively. Consequently, when the individual irreversible reaction Thiele

moduli get near four, g2 can easily be greater than 20. (An exception is Case VII

kinetics, where g2 is identically 0). The other invariant, g3, is frequently smaller in

magnitude than g2 for weakly to moderately diffusion limited reversible reactions, and

tends to be negative for the conventional second-order reactions considered here. When

one or more irreversible reaction Thiele moduli for a problem is near three, either g2

and/or g3 could be larger than twenty in magnitude and convergence of p will take more

than just a few terms. It is best to determine a few of the coefficients in eq. (2.55)

before assessing the unwieldiness of the exact solution. Larger invariant magnitudes

were found to be indicative of steeper concentration gradients and stronger diffusional

limitations. The practical point being made is that for small invariant magnitudes, p is


easily handled with a pocket calculator. A computer evaluation of a lengthy truncated

series is only required for larger invariant magnitudes.

The dependence of p on c"(0) can be very weak for second-order reversible

reactions, since g3 is proportional to powers of the two fundamental parameters added

to similar terms multiplied times c"(0), which is always bounded by one-half for the

reactions being discussed in this chapter. Furthermore, I g3 tends to be a maximum at

the equilibrium conversion, i.e. a root of -Ic"2 + 4 sc" 1 or c'2 + t2c' + '.

For fixed p(1;g2,g3), where g3 is an estimate of the true g3, eq. (2.41) becomes a simple

quadratic equation for c"(0)2

4Wc"(0)2 (50st +12p(1; g2 g3))c"(0) + 6 = 0. 2.58

The situation above is simpler for the two irreversible second-order reactions,

Cases VI and VII, but not better. These reactions should have only two or one modulus

respectively for even cruder analyses. For Table 2.1, Case VI, (A + B -), the two

fundamental parameters depend on only two species' parameters, A,+B and 4+A, as do

the invariants

2 2
B+A A+B 2.59

2 2 3 32 02. 3
( .+B BA)3 +B 2 A+BCA (0) 2 2 2 2.60
g = -- +( B'A-4-B)cA (0)
216 6 3


since CA(0) is also a function of +B, and CI+A. Barring a major difference in

diffusivities of A and B, or a large departure from a stoichiometric surface concentration,

g2 and g3 will be manageable numbers up to quite large P2. For Table 2.1, Case VII,

(2A -- products), there should be only a single parameter, and 4~ and are directly

related as expected. The first invariant, g2 is identically zero and

6 3
S40CA (0) 2.61
g3 27

When g2 = 0, only every third term in the expansion of p is nonzero. Nevertheless, in

both Cases VI and VII there is a limiting computational aspect at large values of the

irreversible reaction Thiele moduli.

This situation is closely related to one degeneracy or simplification of the Pe-

function that occurs when the discriminant of the canonical polynomial, 4o3 g2a g3,

is zero, i.e. g3 27g2 = 0. This is the traditional criterion for a cubic polynomial to

have a repeated root adapted to the canonical form used with the Weierstrass elliptic

function. The discriminant of the canonical polynomial characterizes its roots. When

the discriminant is positive there are three real roots, and when it is negative there is

only one. It is easiest to use the c' equations for further analysis of zero valued


82 2.62
x 27


which can be expanded to

4c'(0)3 + 6 ,,c'(0)2+ 12'Yc'(0) +6~ ,Y 0fst =( -(, 4TY)I- 2.63

This expression can be treated as two cubic equations in c'(0). All six roots have been

found and are given by

o +,, J st- -4 Y 2.64
c'(o) =

-Olst2 Qj 4T 2.65
c'(0) = 2 -4 2.65

where the two roots given by eq. (2.64) are both double roots of eq. (2.63) in addition

to also being the roots of the diffusion modified reaction equilibrium expression, eq.

(2.50). Thus the discriminant equality is expected to be met asymptotically in the limit

as the midplane concentration approaches the diffusion modified reaction equilibrium

concentration. The negative root in eq. (2.66) gives a concentration further from the

surface concentration than the negative root in eq. (2.65), and is physically unachievable.

It is potentially possible that the surface given by the positive root of eq. (2.66) intersects

the c'(0) of (4+,t, 'I) diffusion reaction solution surface, eq. (2.37). This case was left

unexplored, but methods were developed to deal with zero or near zero discriminants in

general in the computer codes used to generate general solution surfaces.

For Table 2.1, Case VII, cA(0) goes to 0 as VA goes to infinity, and g3 goes to

zero, so both g2 and g3 are zero, and again the discriminant approaches zero. When


there is a repeated root, eqs. (2.18), (2.21) and (2.40) can be integrated using

conventional methods to give concentration profiles as functions of either sin2(z) or

sinh2(z) depending on the sign of the repeated root; see below.

One special case in the analysis of the discriminaant occurs when c'(0) = c'

= c' = 0, i.e. a gradient free case, and therefore the second derivative is zero. This

leads to a line in the (4,,, ') parameter space, = 0. Crossing this line in parameter

space corresponds to changing the number of real roots of the canonical polynomial. On

this line, the numerator of eq. (2.37) for c'(z) is zero, and no solutions are obtained for

which c'(0) 0. If the discriminant of the canonical cubic polynomial is zero, and a

double root, ri, exists, then the third root, r3, is equal to -2i,; see Abramowitz and

Stegun (1964). For a positive double real root, ri, p is given by

3 ,
p(z; g2, g3) = fl + > 0 2.66

Since sinh2(x) = (cosh(2x) 1)/2, this gives approximately a first-order concentration

profile expression when sutstituted into either eq. (2.37) or (2.41). The result is that the

dimensionless concentration depends on the inverse of an inverse of a hyperbolic cosine

function plus supporting constants. Two other cases arise

p(z; g2, 3) = P1 < 0 2.67
sin( z)

(z; 2, g)=1 = 0. 2.68

The last equation is for a triple root, which only occurs when all three roots are equal

to zero. These equations give a little more insight into the nature of Weierstrass Pe-

function. The physical significance of the double root cases to the nature of the Pe-

function is that one of the two characteristic periods of p becomes infinite. The Pe-

function then becomes simply periodic. The periodicity of the double poles of the

Weierstrass Pe-function then goes as the inverse of the square of the simple zeroes of

sin(2cx) or sinh(2oiy), for z = x + iy.

If g2 and g3 are both real and if the discriminant is positive, then there are three

distinct real roots, r, > r2 > i3, and p is related to Jacobi elliptic functions, e.g. sn(z),

by relations such as

P1 fs 2.69
P(z; g2, g3) = 3 + 2.69
sn2( r3 z)

Other relations with other Jacobi elliptic functions apply for other root cases. Since no

roots are known in advance due to their dependence on the unknown midplane

concentration, these expressions are both awkward in applications and difficult to

generalize much as was the case for Thiele and his elliptic integrals (1937). Eq. (2.37)

with p is preferred as being more general, since it holds throughout the parameter space

regardless of the number of real roots.


Thus an analytical solution has been derived for the Dirichlet problem for one-

dimensional Cartesian coordinates in terms of the Weierstrass elliptic Pe-function for all

seven cases of Table 2.1, subject only to some practical computational limits as described

above. This solution becomes the focal point for further analysis.

Concentration Profiles From the Analytical Solution

It was mentioned earlier that, if the two invariants both exceed 20, then the

computational effort for determining p becomes significant. Simplifications should occur

in the other extreme of small invariant magnitudes. Consider the problem of determining

c'(0) using eq. (2.37). From eq. (2.56) the expansion of P has the following form

2 4 2_6 _
1 82 g 9 gzz6 3g+g3Z 2.70
(z; g, g) = + ++ + 3g3z + 2.70
P(z; 20 28 1200 6160

In the seven discussed cases, g2 is ab initio known. Therefore the first, second, and

fourth terms in p are known functions of z, while the third and fifth terms are linear in

the unknown invariant, g3. Recall that

4c'(0)3 + 6 t c'(0)2 + 12Yc'(0) + 6~ 'Y 2.71
g3 -216

To find the maximum feasible range of g3, set dg3/dc'(0) equal to zero, giving (after

multiplying by 18)


c'(0)2 + 2,Sc'(0) + T = 0


which is nothing more than the dimensionless diffusion modified rate expression

evaluated at the midplane concentration. So a maximum or minimum in g3 occurs at the

equilibrium concentration. Thus, c'(0), bounded by the surface concentration and

equilibrium concentration, sets the possible range of g3 for given parameters.

Table 2.6 indicates that a five term expansion has 99.9% accuracy at z = 1 when

1g21 < 3 and Ig31 < 3. The five term expansion for p is linear in g3, so eq. (2.37) is

approximately a quartic equation for c'(0). Using just the first two terms of the

expansion for p is within about 10% for this condition (and gets progressively better for

smaller invariants), but leads at once to an easily solvable quadratic equation for c'(0)

-5 0s,, 12 -0.6g2
c'(0) =


where the positive sign must generally be chosen. This equation estimates c'(0) better

than the two term expansion estimates p, and can itself be sufficiently accurate for small

Sg21 and Ig3 |. A two-step method to estimate c'(0) is suggested for I g2, I g3 < 10.

Use the above quadratic equation to estimate c'(0), use the c'(0) obtained to estimate

g3, eq. (2.72), use the five term expansion to estimate P and use eq. (2.37), which

becomes quadratic with an estimate of p, to obtain a new solution for c'(0).

The following example illustrates the procedure. Given

2 1
lst =15 T = 27


calculate g2 from eq. (2.34).

g2 = 9.75 2.75

Find the equilibrium concentration from eq. (2.73)

cq = -2.09167 2.76

For the first step, make a crude estimate of p using

Sz 1 + = 1.49 2.77

Find a first estimate of c'(0) from eq. (2.74) above.

c'(0) -1.21 2.78

Use c'(0) to estimate g3.

g3 = -5.61165 2.79

For the second step estimate p with the five term expansion, eq. (2.71).

(1; g2, 3) = 1.3397 2.80

Recalculate c'(0) with eq. (2.59) at z= 1 and the above p giving

c'(0) -1.945 2.81

Repeating the second step gives c'(0) = -1.948, while the exact solution is -1.946, so

the above estimate is within one per cent. The magnitude of g2 is approximately ten, so

this is a borderline case of the method. The center concentration is at about 93%

potential conversion relative to 100% at diffusion modified reaction equilibrium, so this

example is not in the diffusion free region.


The two-step method is not limited to c'(0) but could be applied to c"(0) or

cA(0). As the magnitudes of the two invariants get larger, the simple hand calculation

above is less accurate. The two-step method for quick hand calculation is thus limited

to the region of small to moderate diffusional effects, even though that may cover 90%

of the possible range of c'(0) values, as in the above example. In the region of strong

diffusional effects, a longer expansion of the Weierstrass Pe- function is required.

First-Order Subcases of the General Problem

For the cases described earlier, there are numerous subcases that exhibit pseudo-

first-order or even true first-order kinetics. Maymo and Cunningham (1966) observed

that for Table 2.1, Case VI, irreversible kinetics, as DB goes to infinity, c goes to cB,,

and the rate becomes approximately equal to the pseudo-first-order form (k ,B A.

Similar arguments, applied to one or more diffusivities, could be extended to any of

Cases I-V.

Cases I-III in Table 2.1 all have another property in that at quite ordinary

parameter values, all the coefficients of the c"2 term, i.e. 'k, vanish as do the a1i in all

of the species' governing ODE's. This arises from the left-hand difference in the

equation for 'i, eq. (2.46). This degenerate first-order behavior occurs not as a limiting

behavior as in Case VI above, but in the middle of the parameter space. This

phenomenon is independent of the initial coordinate geometry (slab, sphere, cylinder),

and is not a statement of reaction equilibrium. This is a degenerate zero in the

discriminant of the canonical cubic polynomial not clearly shown in the development


using c'(0). To explore this situation requires the initial dimensionless species

concentration equations in terms of the aj,. Developing the discriminant in the ai, for a

species, I, gives

g2-27g3 = a [(36a2a3 64ala3) + (108ala2a 54a)a4 (27af)a4 2.82


a4 = -4a,1c,(0) 6a2c(0) -4a3c,(0) 2.83

As seen from eq. (2.82), ai, = 0 also corresponds to a case where the discriminant of

the canonical polynomial equals zero independent of the six c'(0) zeroes derived earlier.

A degeneracy in the definition of c' or c" occurs when a,, = 0. Thus all solutions in

the neighborhood of the plane ai = 0, or I = 0, are strongly first-order in character.

This plane does not apparently correspond to any further unusual behavior of the solution

other than that numerically (a2) = (a2)j, (i,j =A,B,P,R) and 6a2i = D' ,. Something of

the physical significance of 6a2, is thus revealed. In the general second-order reversible

kinetics for Cases I-III, 6a2, can be rewritten as

6a2i = s 12al. 2.84

So for ai, i 0, 6a,, carries the species dependent deviation from first-order behavior for

the coefficients of both the first-order and second-order terms, ci and c?, in eq. (2.18).

The above relationship, eq. (2.84), holds for all species in all seven cases, although its

existence cannot be motivated by the above line of reasoning (search for first-order

behavior), especially for the two irreversible cases where some balance between reactant

and product Thiele moduli is impossible. The initial recognition of eq. (2.84) occurred

prior to the discovery of the two-parameter basis and was the first clue that more

fundamental parameters might exist. When ali, or I, is zero, it makes little sense to use

c', defined with muliplication by ali, or c", defined with division by 'I. For Cases I-

III, when a,, = 0

d2ci 2.85
2 = 6 a2 c, + 2a3, = Q c + 2a3s 2.85

which has the solution

2 12a
c,(z) = s + 2as cosh(Osz) 2 2.86
2 Csth(4 ) o s t 2
1st ( i st( 1st

with the center concentration given by

ci(0) = 0st + [2a3 ( ist s os 2.87

This form of the expression permits ready comparisons with the similar expressions for

cylindrical and spherical coordinates, Table 2.7 and 2.8. (The functions 4 2cosh(4,t),

sinh(c(P,), and tlo(4L) in Table 2.7 all have series expansions in odd powers of 4,

only, beginning with the same first-order term.) In the limit of small V,, a3i also goes

to zero, thus c,(O) goes to one in this limit as expected. For some neighborhood of a1i

about zero, this expression can be used to estimate ci(0) for second-order cases. Tables

2.7 and 2.8 also give the dimensionless concentration profiles and center concentrations

for the Robin's problem in all three coordinate systems. The general Robin's problem

is discussed in Chapter Four, but solutions for the first-order subcase are included here

for convenience. The internal effectiveness factor is also defined more formally than in

Chapter One by

(p + 1)D.d., 2.88
lint = L2( is, js'" )

where p is zero for a slab, one for a cylinder and two for a sphere; see Aris (1975) or

Lee (1985) for general details. For the slab-like medium, the solution is

tanh( Q )
l int anh( a, =0, T' =0 2.89

which is identical in form to those for simple first-order irreversible or reversible kinetics

in Cartesian coordinates, Chapter One. When ali = 0, an analytical solution to the

reduced first-order kinetics problem can be derived in spherical or cylindrical coordinate

systems as well. Internal effectiveness factors for the three common geometries are

given in Table 2.9 for the Dirichlet and Robin's problem.

It is worth noting that, given fixed surface concentrations, temperature could, in

theory, be manipulated for Cases I-III, such that a1i = 0, through temperature's effect

Table 2.7 Degenerate First-Order Concentration Profiles

Dirichlet Problem ci.

Slab ~ S + 2a3. 2a
coh(i ) 3J cosh(40 1z) 2
0cs(cosh(l i 05

Cylinder 0 + 2a3, 2a3,
2 0(01z) 2
Io- lst ) 2 "

Sphere 'Lst + 2a3, sinh($tz) 2a3
0 sinh(4D) 0 z 2
Robin's Problem

Bim(l + 2a3i/ O.,,)
Slab 1i,,sinh( O,) + Bicosh(i,,,)

-cosh(01,,z) 2a3,/ st

Bim(1 +2a3t/ Ct)
Cylinder 1st I (01st) + Bimlo( 0 ,)

*Io(1$,tz) 2a3i,/L 0

Bim,(1 + 2asI 2
(Bim 1)
cosh(Qlst) sitih(N/ st )
Sphere cosh( + st ih

sinh(4i1,z) 2a3,
1stZ 1st


Table 2.8 Degenerate First-Order Midplane Concentrations

Dirichlet Problem i

Slab lt + 2a3i[O s, cosh( )l1l / 02st
Slab 1sst Is

Cylinder lr2 isto (air s
1stIo ( Ist,)

Sphere 0t + 2a3,[ s ch(s)] / ist

Robin's Problem

1 +2a3il/ st 2a3,

Slab (+s2 1 2a
1 + 2a3i/01lst 2a3i
Cylinder Io(1) I ( t) 'Lit

2 2a3,
Shr1 + 2asDis,
sinh(ez,,) cosh(ol,,) sinh(ei,) 4D
Sphere 01st Bim 4 IstBim

on the forward and reverse rate constants, forcing first-order behavior to occur. No

obvious benefit results.

Internal Effectiveness Factors of Second-Order Reactions

The internal effectiveness factor measures the ratio of actual consumption of a

species inside a medium relative to the consumption possible in the absence of

concentration gradients, i.e. at surface compositions. For the single reaction systems

considered in this article, this expression is given in numerous texts including Aris

(1975) and takes the following form for second-order systems in terms of c' or c"

(dc' (dc"\
dz )l dz= 2.90
'lint m 1


(dc' _2c,(0)3 2 2 2.91
= 3- 2c'(O ,c'(O) 2c'(0) 2.91
Sdz )i s 3

where the sign is chosen to match the sign of I, and where

( dc" 2'Pc"(0)' 3 c( 2 c0) 2.92
dc 2 0 1c "(0)2 + 2c(0) 2.92
dz i \ 3

Table 2.9 Degenerate First-Order Effectiveness Factors

Dirichlet Problem




Robin's Problem





lstlo(Ou )
3 3



1Qi13 13)



3coth(d,,,) 3/1st

+ 0Stcoth(-lt) 1

Alternatively, the dimensionless species' concentration eq. (2.88) can be used

dci =, /4a1i(1 ci(0)) + 6a2,(l c(0)) + 4a3,(l -c,(O)) 2.93

The remaining portion of eq. (2.88) for the internal effectiveness factor can also be

expressed in terms of the aj. The external dimensionless rate occurs at ci = 1, the

surface concentration, in the right-hand side of eq. (2.18). It can be shown that the

internal effectiveness factor is given by

(dc, )
+1)dz 2.94
= ____ \dz 2
Tlint 6as, + 6a2i + 2a3,

Table 2.10 gives the expression in the denominator for each species in all seven cases

for the slab-like geometry (p = 0). It should be observed that although ali, a2i and a3,

vary from species to species within a case when using these expressions, I,, and do

not, and r can be determined from these alone, rather than the a,j! No species

dependence is suggested for it.

Methods exist for estimating the internal effectiveness factor based on the

generalized Thiele modulus approach, Lee (1985), Aris (1975). The asymptotic form of

the internal effectiveness factor based on a generalized Thiele modulus is often

sufficiently accurate for engineering requirements. This expression matches the

asymptotic behavior at both the diffusion free and diffusion dominated limits, but is

Table 2.10 Denominators for the Internal Effectiveness Factor Expression

Case Species (dci/

I A B +A R+P/rAP
C A p2 C

I P + IAP +A


II A 242A 24p+/rAP

II P p rAP2A

II R P+R A2A +R +P

III A 2.t 24p /rAP

III P 2p rAP A

IV A -+A- ~P rAP

IV B IA+B- 'IA+BP'/ rApB+A


V A 2IA 24/rAP

V P p rA 2A

VI A $+A


VII A 2 2A


uncertain, and not guaranteed to be even approximately correct, in intermediate ranges.

Conventional approaches to estimating midplane concentrations are observed to fail

badly, although limiting asymptotes of the internal effectiveness factor are usually

approximately correct. This is unfortunate, since the internal effectiveness factor for the

slab geometry can be written as an algebraic function of the ab initio known parameters

plus the midplane concentration. The slab midplane concentration is immediately

available from the new analytical solution presented here, plus the analysis suggests

several new ways to estimate it under various limiting situations.

Separate plots of the internal effectiveness factor versus the relevant irreversible

reaction Thiele moduli can be prepared for both Case VI and Case VII kinetics. Because

these two irreversible reaction systems have been studied by others numerically, these

plots will be presented here using the analytical solution. The results for these two cases

are also part of the general plots of the internal effectiveness factor versus l,, and *I

presented later. Internal effectiveness factors have been obtained numerically by Bischoff

(1965) for the problem

Dd k'cn 2.95

where n is the reaction order and k' is a rate constant. This plot is given in many

standard heterogeneous kinetics texts. The comparable problem in this paper, Case VII

of Table 2.1, is

Dd2 = 2k,f2 2.96
A dx2 CkA

with the two fundamental moduli given by

4L 4L 4 bk2 A2
s 4Lf As 2.97
t D 4 2o

There is really only a single parameter, since 41,t and P are related. It is almost

certainly the best known simultaneous second-order reaction and diffusion problem result.

Plots of qit are made in Bischoff versus a generalized Thiele modulus, $o, given by

2 L2(n+1)k' ,- L23(2kf)A 2.98
2G = D C 2DA45A 2.981A
20 2,A

modified here to match the current nomenclature and to coincide with second-order

reaction systems. A feature of the generalized Thiele modulus development is that

1 int I o 2.99

Noting that k' = 2kf, substitution of eqs. (2.97) and (2.98) into eq. (2.99) yields

Sint Ol r 00 2.100
V3 st

as the equivalent expression. Fig. 2.6 compares the exact solution for the irreversible

power law kinetics of Case VII with an asymptotically correct expression based on the

generalized Thiele modulus approach. The asymptotic curve is based on

9 4-












'4 ^