M. J. SHAH
BURTON DAVIDSON
KINE TIC MODELS FOR CONSECUTIVE HETEROGENEOUS REACTIONS The equation is semitheoretical-when only one of the steps considered is truly ratein the field of heterogeneous controlling, it degenerates catalysis beget sophisticated to the exact equation which o p t i m i z a t i o n problems. applies. This interpolaThere remains a basic need tion method provides, in a for realistic, useful kinetic rational manner, an accurate rate equations, and therefore mathematical simulation of a need for methods of dethe complex chemical reacveloping such expressions tions taking place. On the from experimental data. I n other hand, the number of heterogeneous catalysis, reparameters in nonlinear action can take place in form is not changed, and the several consecutive steps, or difficult problems involved in parallel steps. Reaction in nonlinear searches are rates may therefore be reptherefore not avoided. Furresented by equations of thermore, the exact physical tremendous complexity, and meaning of the interpolasolution of such equations tion parameters (including s) would be at least a formidaPRESSURE cannot be easily interpreted. ble task. These liabilities reflect the One way out of this mathempirical side of the method. Any change in process ematical maze is to assume a simple physico-chemical conditions which would affect the type of reaction mechmechanism-to assume that a single step controls reaction anism prevailing would require only minor adjustments rate. This route has been followed in the past, of in the interpolation model structure. These changes necessity, and results in rate equations such as Equation would be, for example, assigning the value of zero to cer1. But the concept has probably been overworked as tain r.)s in the equation. far as heterogeneous catalysis is concerned, and is, in fact, a carry-over from the precomputer era. Reaction Mechanism Consideration of more than one rate-controlling step, Model selection has traditionally been based on Le., consecutive heterogeneous reaction, leads to cumberLangmuir-Hinshelwood activated adsorption-surface resome rate expressions. A mathematical model for action mechanisms, which have been popularized by a relatively simple case with two rate-controlling steps, Hougen and Watson (8). According to their theory, a as dweloped by Bischoff and Froment (I), is given as series of specific events must take place: Equation 2. Evaluation of the nonlinear parameters of -physical mass transfer of reactants and products to such a model is impractical for most problems. and from the gross exterior surface of the catalyst parOne is therefore led to another route-the development ticle and the main body of the fluid interpolation expressions which require fewer difficult --diffusional and flow transfer of reactants and prodalgebraic manipulations. For consecutive heterogeneous ucts in and out of the pore structure of the catalyst reactions, we suggest an expression of the form: particle where reaction takes place at interior interfaces r, = C (r.a))''s -chemisorption of reactants at active sites --surface reaction of adsorbed reactants to form The effect of such a method on fitting a curve to experichemically adsorbed products mental data, is illustrated above. A compromise is found -activated chemidesorption of products at the catalyst between the two extremes represented by the simple interface models. procedures for Design complex reactor systems
18
EXPERIMENTAL DATA MAY FALL INTO A TWILIGHT ZONE, NOT DEFINED BY ANY SINGLE RATE-CONTROLLING STEP
INDUSTRIAL AND ENGINEERING CHEMISTRY
Interpretation of experimental data in terms of this basic mechanism, involves three steps : model selection, evaluation of parameters in the selected model, and simplification of the model expression. T o illustrate the various techniques used, we will use the experimental data given in Table I throughout this article. Model Selection
Perona and Thodos ( 7 4 , Thaller and Thodos (20), and Ford and Perlmutter (7) have studied the brasscatalyzed reaction involving the dehydrogenation of sec-butyl alcohol to methyl ethyl ketone and hydrogen from the traditional Hougen-Watson viewpoint. The various rate equations used by these authors can be found elsewhere (7, 22). The selection of the various models which led to these expressions was based on the same premise-it was assumed that only one type of chemical dissociation mechanism and only a single rate controlling step were involved in any one of the reaction rate models chosen. A typical example is Equation 1 . For the derivation of this equation, the reaction was assumed to take place in three steps : (1) adsorption of the alcohol on a single active site ( 2 ) surface reaction between the adsorbed alcohol and an adjacent active site to produce the adsorbed products (3) desorption of the products The first and third steps were assumed to be at equilibrium, with the second step controlling the rate of the reaction. The partial pressures were written as interfacial partial pressures which can be estimated using conventional methods (6, 8, 22). This essentially accounts for physical mass transfer effects. From the
Yang and Hougen tables (ZZ), Equation 1 can be constructed without precipitating tedious algebraic steps. Until recently, nearly every chemical kinetic paper dealing with the Hougen-Watson approach selected reaction rate models similar to Equation 1 . More involved models may be developed, however, by assuming that more than one step is rate-determining, while the type of surface-chemical mechanism remains the same. This is commonly referred to as consecutive heterogeneous reaction. For example, Bischoff and Froment (7) have developed a n analytical kinetic model for the dehydrogenation of sec-butyl alcohol reaction based on two rate-determining steps, which are: the surface reaction between chemisorbed alcohol and an active site, and the desorption of chemisorbed hydrogen. Equation 2 is the resulting rate equation. Use of Initial Rates
For the purposes of mechanistic studies, Equations 1 and 2 are so complicated that the evaluation of the various nonlinear rate and adsorption equilibrium constants from the experimental data, even by means of nonlinear least-squares estimation, is difficult. I t is therefore natural to use initial rate data to simplify the anticipated equations. This procedure eventually eliminates the products from the variables in the models. Equation 1 would then become Equation 3. If desorption of hydrogen is controlling with a dualsite surface reaction mechanism, then the initial rate expression is given in Equation 5, which is obtained from Equation 4 by substitution of p K = p H = 0. If Equation 2 is simplified in a similar manner, Bischoff and Froment (7) showed that the initial rate expression,
SINGLE RATE-CONTROLLING STEP
TWO RATE-CONTROLLING STEPS
SURFACE REACTION CONTROLLING:
SURFACE REACTION AND PRODUCT DESORPTION CONTROLLING:
INITIAL RATE EQUATION IS:
TAi
=
L J ~ R ~ A ~ A
(l f
KAPA)'
DESORPTION O F HYDROGEN CONTROLLING: kHKL
rA =
P K (1
[@A
-
(PK$H/K)
KAPA 3- KKPK
+
c'
1
K K H
-
=
kH2K
CORRESPONDING INITIAL RATE EQUATION IS:
INITIAL RATE EQUATION IS: TAi
-
kHL/KH
V O L . 57
NO. 1 0 O C T O B E R 1 9 6 5
19
AN UNPROFITABLE APPROACH INTERPOLATION F O R M U LA
ALTERNATIVE FORMULATION OF EQUATION 6, F O R T H E EXAMPLE GIVEN:
SHAH-SYN SEARCH TECHNIQUE APPLIED T O EQUATION 12
min of
(rAi
- $),"
On2
=
n
n
where D, = (rAi - q5)n, s = -0.7 GENERAL INTERPOLATION FORMULA
when both desorption and surface rates are controlling with a dual-site mechanism, reduces to Equation 6. Furthermore, it was shown that when the desorption of hydrogen is controlling ( k , + a), the initial rate expression from Equation 6 would become: 'Ai
=
kH
(7)
This contains the factor L/'KH. If, on the other hand, kH +- m , surface reaction is controlling, and Equation 3 is obtained with the factor Ls implicitly contained within the kinetic constant, k,. Other more complicated models than Equations 2 or 3 can be derived by considering other rate-controlling steps and other chemical reaction dissociation steps, For example, a single model can be selected based on the three mutually controlling steps of surface reaction (dual-site), desorption of hydrogen (dual-site) in parallel with surface reaction (single-site), and vapor phase dissociation. This model would indeed yield a reaction rate equation of tremendous complexity; and fitting experimental kinetic data to such an expression would certainly not be easily accomplished. Evaluation of Parameters in Selected Reaction Models
The estimation of rate constants and adsorption equilibrium constants in kinetic models has been a subject of speculation because of the nonlinear nature of the derived rate equations. Before the age of the high speed digital computers, the parameters in nonlinear models, similar to Equation 1 , were evaluated using linearized techniques described by Hougen and Watson (8). Linearization permitted the use of either graphical or ordinary least-square techniques. Numerous examples of these methods may be found in the literature (6, 7, 74, 20). Chou ( 4 ) has revi,ewed the limitations of the least-squares method and Opfell and Sage (73) have discussed a variety of applications of least-squares procedures to other types of equations with both linear 20
INDUSTRIAL AND ENGINEERING CHEMISTRY
and nonlinear parameters. Blackmore and Hoerl (2) have thoroughly treated the applic,ation of least squares to nonlinear models of the type given by Equation 1. They have demonstrated the pitfalls associated with forced linearization of nonlinear models and have given quantitative evidence of the effects. Although there TABLE
I. I N I T I A L REACTION RATES-DEHYDROGENATION O F sec-BUTYL ALCOHOL
Feed: 100% sec-butyl alcohol Catalyst: 1 gram of brass spheres. (50- to 60-mesh) mixed with glass spheres of the same size in the proportion of 20 parts of glass to 1 of catalyst
Tpm$, ,
Pressure,
F.
Atm.
Feed Rate, Lb. MoleslHr.
1 .o 7.0 4.0 10.0 14.6 5.5 8.5 3.0 0.22 1 .o 1 .o 3.0 1 .o 3.0 5.0 7.0 9.6 2.0
0.01359 0.01366 0.01394 0.01367 0.01398 0.01389 0.01384 0.01392 0.01362 0.01390 0.01396 0.01392 0.01411 0.01400 0.01401 0.01374 0.01342 0.01386
O
600 600 600 600 600 600 600 600 600 600 550 550 575 575 575 575 575 550
PA,
r A i , Initial Rate, Lb. Moles Alcohol/Hr.Lb. Catalyst
0.0392 0.0416 0.0416 0.0326 0.0247 0.0415 0.0376 0.0420 0.0295 0.0410 0.0115 0.0161 0.0227 0.0277 0.0255 0.0217 0.0183 0.0146
6 8 IO 12 14 16 18 PRESSURE, ATMOSfHERES
0 2 4 F i p e 7 . Pbr of Equatim 6 at 575' and 6oooF. w'fh pormnekrs cddated from Equafion8
F i p t 2. Plot of Equafion 12 at 575O an SJO' F. u ith parmfers cdculoted from Eguohon 73, and s = -0.7
does not seem to be any general theory for estimating parameters in nonlinear models using the least-squares criterion on the nonlinear model (72), Rubin (76) and others (2, 70) have outlined various useful computat i 0 ~ 1schemes for fitting nonlinear reaction rate equations to data. Some recent methods have been described by Peterson (73, Lapidus and Peterson (77), Shah and Syn (79), and Davidson (5) for nonlinear examinationof rate constants. SimplillcoHon of R a t . Expressions
As stated earlier, Equation 1 represents one method of simplification. A further simplification, popularized by Weller (27), is based on a purely empirical powerlaw formulation that is designed for expediency for we in design calculations. T h e apparent pitfalls of this method have been discussed by Davidson and Thodos (6),Boudart (3),and others ( 5 2 2 ) . Several other methods of solving for the parameters of a model with two rateeontrolling steps have been tried. For illustration, the original experimental data of Thaller and Thodos (20) for initial reaction rates of dehydrogenation of sec-butyl alcohol have been used. The data are shown in Table I.
M. J . S h h is with the Data Rocessing Section, General Products Diuision of International Busincss Machzncs Cap., San Jose, Calif. Burton Davidson is Professor of Chemical Engiwering, at Rutgns, the State University, New BnuuWuk, N. J . The assistance of W. M . Syn with nwnnual 'calculationrisgnrrcfurry acknowledged. Burton Davidsm wishes to atknowledge thejnancial assistance of the Research Council of Rutgns, and a San Jose Collegefaculty grant. AUTHOR
Figwc 3. Plof of Equaf~~ns 6 and 12 at 6oooand 57.5' F.
The model selected was Equation 6. Estimation of the nonlinear parameters (kR, k,, and KA) was a p proached initially from the traditional least-squares method, following Equation 8 :
The corresponding normal equations were obtained by setting the iirst partial differentials with respect to k,, k,, and K A equal to zero. The first partials are sufficient for the analytical determination of the parameters, but the second partials have to satisfy Lagrange's criteria to guarantee a minimum. The resulting normal equations, Equations 9, 10, and 11, were derived by Hughot and Davidson (9). An analytical solution of Equations 9 to 11 is infeasible in closed form. A numerical solution was attempted by Hughot and Davidwn (9) on an IBM 1620 using a Fortran Search program, but the results, unfortunately, were sketchy and inconclusive. A second attempt was tried using the d m t search technique of Shah and Syn (79) on Equation 8. Because one is forced to employ a numerical procedure for estimation of the rate constants in the normal equations, it is considerably less cumbersome to employ Equation 8 rather than Equations 9 to 11. The search technique of Shah and Syn is suitable for evaluation of parameters in nonlinear formulations such as Equation 8. Often in nonlinear estimation problems of this nature, determination of the minimum squared error leads to more than one minimum (polymodality), so that the numerical search may terminate with erroneous estimates of the parametus. The technique of Shah and Syn was programmed in Fortran on an IBM 7090 and i s particularly VOL 5 7
NO. 1 0 O C T O B E R 1 9 6 5
21
well suited to polymodal search. As far as the technique is concerned, up to four parameters may be estimated simultaneously with alacrity and the degree of nonlinearity- or complexity is immaterial. For more than four parameters, as in Equation 2, special sequential techniques have to be employed if the Shah-Syn (19) search technique is to be used. The data of 575' and 600' F. from Table I were used in the nonlinear least-squares problem given by Equation 8, and were processed to yield estimates of k,, kH, and K A . The data at 550' F. were considered inadequate for least-squares analysis and were not used. The results are shown in Tables I1 and 111. Theoretical reaction rates as predicted by Equation 6 with the values of k,, k,, and K A from Equation 8 are plotted us. pressure in Figure 1. The actual experimental data from Table I are included for comparison. The agreement between calculated and experimental rates is excellent. Improved correlation at the higher pressure end of the spectrum may be anticipated from consideration of additional rate-controlling steps ; and parallel as well as consecutive mechanisms may be significant. All parameters are positive as required by theory. The possibility of further simplification of the foregoing illustration is apparent in view of the complexity of the rate expressions for heterogeneous kinetics when more than one rate-controlling step is taken into consideration. The determination of rate constants from these expressions becomes even more laborious ; and as the number of these constants increases, it becomes harder to extract meaningful values of the constants from the data. Lapidus and Peterson ( 7 7 ) have described techniques for estimation of the constants by trying out various rate expressions, each derived for a single ratecontrolling step. T h e expression that will give the best fit may perhaps be taken as the satisfactory one. However, Bischoff and Froment ( 7 ) have pointed out how this may lead to erroneous mechanisms. Interpolation Method
There still remains a need for a technique, for design purposes as well as for a quick determination of rate controlling mechanisms, that can predict accurate rates and realistic estimates of the kinetic constants from experimental data. An interpolation formula is proposed which not only simplifies computational problems, but also takes into consideration several rate-controlling mechanisms in a natural manner. The method is quite accurate for extrapolation in design calculations. I n the illustration of the brass-catalyzed dehydrogenation of sec-butyl alcohol, the expression for initial reaction rate, Equation 6, simplies to Equation 7 if hydrogen desorption is rate controlling ( k R + m ) and to Equation 3, if surface reaction rates are controlling (kH + m ) . An interpolation formula when both rate-controlling steps are equally important is expressed in Equation 12. Equation 1 2 is formulated as an alternative to Equation 6. The constant s need not be a n integer. Now if s is negative, Equation 1 2 will reduce to Equation 7 when k, + m and to Equation 3 when k H + a . Thus, the 22
I N D U S T R I A L A N D E N G I N E E R I N G CHEMISTRY
interpolation formula will still be valid when only a single rate-controlling condition exists. T o indicate the accuracy of Equation 1 2 relative to Equation 6, the nonlinear least-squares problem shown in Equation 13 was formulated. The data from Table I were used and the nonlinear parameters ( k H , k,, and K A ) were estimated according to the Shah-Syn technique. The results are tabulated in Tables I V and V. A plot of Equation 1 2 is given in Figure 2 to indicate agreement of experimental and calculated values. The parameters were determined from Equation 13. The close agreement between Equation 6 and Equation 1 2 relative to the experimental data is shown in Figure 3. This fact, coupled with the close agreement between calculated parameters from the different equations, suggests grounds for a generalization of Equation 12, as formulated in Equation 14. I n Equation 14, the subscript n indicates a reaction kinetic expression involving only a single controlling step of the type given by Equations 1 and 5. The exact form of the rn functions may be determined, without derivation, from the Yang and Hougen tables (22). If initial rate is needed, then rn represents an initial TABLE II.
NONLINEAR LEAST-SQUARES PARAMETER E S T I M A T I O N F R O M EQUATION 8
kx X 702, Lb. Moles T,
O
F.
_ _
Alcohol/Hr.Lb. Catalyst
k~ X 702, Lb. Moles Alcohol/Hr .Lb. Catalyst
K A X 702, Atm.-'
575 600 WI
H PARAMETERS F R O M EQUATION 8
Rates, Lb. Moles Alcohol/Hr.Lb . Catalyst Temp., F.
Pres?Aa X 1 (rat, sure. 7.4ia) Atm. Measured ! Calcd. x 102 _ _ _ _ _ _ _ ~
1 575
600
5 7 9.6
0.22 1 .oo 3.00 4.00 5.50 7.00 8.50 10.00 14.60
~
2.27 2.55 2.17 1.83
2.95 3.92 4.20 4.16 4.15 4.16 3.76 3.26 2.47 1
2.40 2.58 2.41 2.22 2.01
4.10 3.92 3.75 3.59 3.44 3.07
~
1
Dedation,
70
-0.13 0.19 0.14 - 0.05 - 0.07
5.73 6.86 5.49 2.30 9.29
Average
5.94
0.06 0.23 0.4 1 0.17 -0.18 -0.60
~
Average
1.02 6.13 0.23 1.44 5.54 9.86 4.52 5.52 24.30 6.51
reaction rate expression similar to Equations 3 and 7. As far as design procedures are concerned, Equation 14 is not only simpler than the expressions of the type given by Bischoff and Froment ( 7 ) but it also takes into consideration several consecutive rate-controlling steps in a natural manner (carte blanche). The only condition in Equation 14 is that s must be chosen so that, for single rate-controlling conditions, the expression must degenerate to the rate expression for this condition. A numerical value for Is/ of ”8 to ”4 is recommended based on prior successful experience (79) with these kinds of interpolation formulations. I t is best to guess at s prior to parameter estimation. A value of s = -0.7 seems to work well with the HougenWatson rate equations. Conclusions
Until derivations of analytical expressions for more than two rate-controlling steps are available, the validity of Equation 1 4 cannot be fully ascertained. Furthermore, one also needs to consider the utility of the interpolation formula for the case where both parallel and consecutive mechanisms are controlling. TABLE IV.
NONLINEAR LEAST-SQUARES PARAMETER E S T I M A T I O N FROM EQUATION 13
I
Parameter k H x To2, Lb. Moles
OF.
575 600
I
x To2,
NOMENCLATURE
Lb. Moles Alcohol/Hr.Lb. Catalyst
Alcohol /Hr. Lb. Catalyst
Temp., T,
kR
A potentially useful application of the interpolation formula is associated with mechanism screening for systems displaying more than one rate-controlling step. The situation may arise whereby several complicated, dual, rate-controlling mechanisms of the type given by Equation 2 or G would have to be derived and then subjected to nonlinear curve-fitting analysis. The number of such expressions may be reduced by performing the nonlinear curve-fitting analysis on easily derivable interpolation equations, followed by the criterion that all parameters calculated by this procedure must be positive. Models yielding negative constants are then discarded, and the procedure is repeated until at least one reaction rate model produces positive constants. I t would then be necessary to derive the theoretical expression (similar to Equation 2 or 6) for this model followed by a suitable nonlinear least-squares curvefitting analysis. Following the suggestions of Blackmore (2) and others (77)) it would be advisable to evaluate simultaneously the temperature dependence of the parameters in a single iteration routine using the Arrhenius equation and the theory of absolute reaction rate equations in situ with all of the temperature data. This procedure yields the optimum values for activation energies, heats of activation, and entropies of activation, providing, of course, that proper experimental design and accurate data are used (2, ?O,7 7 ) .
K A X TO2, Atm.-l
Ka = adsorption equilibrium constant for sec-butyl alcohol KR K
o::l
kH kE
CALCULATED RATES FROM EQUATION 12 w17 I PARAh ETERS FROM EQUATION 13 Rates, L 6. Moles Alcohol/Hr Lb. Catalyst
TABLE V.
.-
Temp., ’ F.
575
600
Pressure, Atm.
TAi
x
TO2, Measured
?-Ai
x
(rAm.
-
702, TAic) X ?02 Calcd. _______
Deviation,
-0.08 0.1 1 0.08 0.02
3.52 3.99 3.14 0.92
9.6
1.83
1.99
-0.16
8.04
3.00 4.00 5.50 7.00 8.50 10.00 14.60
4.20 4.16 4.15 4.1 6 3.76 3.26 2.47
4.35 4.21 3.97 3.74 3.53 3.35 2.88
Average
3.93
0.27
9.15
- 0.32
0.18 0.42 0.23 -0.09 -0.41
8.1 6 3.57 1.20 4.34 0.10 6.12 2.76 6.60
Average
6.67
-0.15
- 0.05
1
rAi
REFERENCES
2.35 2.66 2.47 2.15
2.68 4.24
rA
s
2.27 2.77 2.55 2.17
2.95 3.92
L
pa p, pK
%
1 3 5 8
0.22 1.00
adsorption equilibrium constant for hydrogen adsorption equilibrium constant for methyl ethyl ketone = surface reaction equilibrium constant = overall thermodynamic equilibrium constant = rate coefficient for hydrogen controlling = rate coefficient for surface reaction controlling = total active sites per mass of catalyst = interfacial partial pressure of sec-butyl alcohol = interfacial partial pressure of hydrogen = interfacial partial pressure of methyl ethyl ketone = rate of reaction based on sec-butyl alcohol per mass of catalyst = initial rate of reaction = interpolation constant or number of adjacent active sites
KH = KK =
(1) Bischoff, K. B., Froment, G. F., IND.ENG.CHEM.FUNDAMENTALS 1 (3), 195-200 (1962). (2) Blakemore, J. W., Hoerl, A. E., Chem. Engr. Sympostum Ser. 59 (42), 14-27 (1963). (3) Boudart,M., A.I.Ch.E.J.2, 62 (1956). (4) Chou, Chan-Hui, IND.ENO. CHEM.50, (51, 799-802 (1958). (5) Davidson, B., “Kinetics of the Catalytic Oxidation of Sulfur Dioxide,” Ph.D* dissertation, June 1963, Northwestern University. (6) Davidson, B., Thodos, G., A.1.Ch.E. J. I O , 568 (1964). (7) Ford, F. E., Perlmutter, D. P., Chem. Eng. Sa.19, 371-378 (1964). (8) Hougen, 0. A,, Watson, K. M., “Chemical Process Principles,” Part 3, Wiley, 1952. (9) Hughot, E., Davidson, B., private communication. (10) Hunter, W. G., Kittrell, J. R., Watson, C. D., Abstract 18d., presented a t the 57th Annual Meeting of A.I.Ch.E., Boston, December 1964. (11) Lapidus, L., Peterson, T. I., Abstract 18c, paper presented a t the 57th Annual A.1.Ch.E. Meeting, Boston, December 1964). (12) Linnik, Y. V., “Method of Least Squares and Principles of the Theory of Observations,” Pergamon Press, New York, 1961. (13) Opfell, J. B., Sage, B. H., IND.ENG.CHEM.50 (5), 803-806 (1958). (14) Perona, J. J., Thodos, G . A., A 1.Ch.E. J. 9, (2), 230-235 (1957). (15) Peterson, T. I., Chem. Eng. Sa‘. 17,203 (1962). (16) Rubin, D. I., Chcm. Engr. Pror. Symp. Ser, 59 (42), 90-94 (1963). (17) Shah, M. J., Ph.D. dissertation, University of California, Berkeley, 1961. (18) Shah, M. J., Tech. Rept. TR-02.273 10-01-63, I.B.M., San Jose, Calif. (19) Shah,M. J., Syn, W. M., Int. J . Camp. Math. (1964). (20) Thaller, L. H., Thodos, G . , A.1.Ch.E. J. 6 (3), 369-73 (1960). (21) Weller, S.,A.1.Ch.E. J. 2, 59 (1956). (22) Yang, K. H., Hougen, 0. A,, Chem. Eng. Progr., 46 (3), 146-157 (1950).
VOL. 5 7 NO. 1 0 O C T O B E R 1 9 6 5
23
