Ind. Eng. Chem. Res. 1991,30,639-645 Registry No. NH3, 7664-41-7;HN03,7697-37-2;furfural, 98-01-1.
Literature Cited Aiken. R. G . Humic substances in soil. sediment and water:. Wilev: New York, 1985. Box. G. E. P. Estadietica Dara inuestigadores: Revert& Barcelona. 1988.
Caro, N.; Frank, A. R. German Patent 559,254,28,1932. Chakrabartty, R. K.; Haque, R. Production of furfural from agroindustrial wastes under normal pressure. Indian J. Technol. 1975, 13 (7),315-318. Coca, J.; et al. Production of nitrogenous humic fertilizer by oxidation-ammoniation of lignite. Znd. Eng. Chem. Prod. Res. Dev. 1984,23 (4),620-624. Coca, J.; et al. Oxiammoniation of pine-bark particles. Can. J. Chem. Eng. 1985,63(5),835-839. Coca, J.; et al. Production of furfural by acid hydrolysis of corn-cobs. Accepted for publication in J. Chem. Technol. Biotechnol., 1990. Grethlein, H. E. Statistical design of experiments for optimizing the casting variables for cellulose acetate membranes. In Reverse osmosis and synthetic membranes theory-technology-engineering; Sourirajan, S., Ed.; NRCC: Ottawa, 1977. Giiroz, K. Oxy-ammoniation of Elbistan lignite to produce a nit-
639
rogenous fertilizer. Fuel 1980,59 (ll),772-776. Hokanson, A. E.;Katzen, R. Chemicals from wood wastee. Chem. Eng. B o g . 1978,74 (l), 67-71. Huibers, D. T. A.; Jones, M. W. Fuela and chemical feedstocks from lignocellulosic bioknass. Can. J. Chem. Eng. 1980,58,71%722. Kim, Y. K.; et al. Fertilizer from the oxidative ammoniation of sawdust. Znd. Eng. Chem. R o d . Res. Dev. 1981,20 (2),205-212. Kizer, 0.; et al. Etude expdrimentale de l'oxydation catalytique du benzene en anhydride maldique on couche fluidisde. Chem. Eng. J. 1977,14,205-215. Martinez, A,; et al. Programas BASIC para la optimacibn de funciones. 3. Mdtodo de Powell. Zng. Quim. 1985,8,97-101. Official Methods of Analysis. Horwitz, William, Ed.; AOAC: Washington, DC, 1975. Popa, V. I. Investigations in the field of lignine valorification. V. Products of oxidative nitration. Cellulose Chem. Technol. 1985, 19,657461. Stevenson, F. J. Humus Chemistry: Genesis, composition, reactions; Wiley: New York, 1982. I Toynbee, D. A.; Fleming, A. K. Air oxidation of subbituminous coal. Fuel 1963,42,379-387.
Received for review January 22, 1990 Revised manuscript received August 27, 1990 Accepted September 17, 1990
Difficulties in the Application of Sequential Experimental Design to Kinetic Model Discrimination and Parameter Estimation L u c h e z a r A. Petrov,* Alexander
E.Eliyas, a n d Christ0 S . Maximov
Institute of Kinetics and Catalysis, Bulgarian Academy of Sciences, str. Acad. G . Bonchev 61-11, Sofia 1113, Bulgaria
Sequential experimental design has been applied to kinetic model discrimination and parameter estimation, based on data on uninhibited C2H40oxidation over a Ag catalyst. The algorithm of Box and Hill and the generalized criterion of Box, Hill, and Wichern have been used. Difficulties in the practical realization of the sequential design of kinetic experiments have been discussed, and a new approach for its application has been proposed. The approach is based on discrimination over a discrete data field, fded up uniformly, after removal of the outliers and preliminary estimation of the parameters, using the complete data set. Introduction A great number of studies have been published in which different variations of the method of sequential experimental design have been proposed for the discrimination between kinetic models and parameter estimation (Agarwal and Brisk, 1985,Box and Hill, 1967;Box and Hunter, 1965; Box and Lucas, 1959;Froment and Hosten, 1981;Hosten and Froment, 1976; Hsiang and Reilly, 1971;Box and Hunter, 1963; Hunter et aL, 1969;Hunter and Reiner, 1965, Kafarov et al., 1974;Reilly, 1970;Reilly and Blau, 1974). However, these methods have been applied to real experimental data only in a few studies (Bibin, 1987;Rozycki, 1987;Dumez and Froment, 1976;Dumez et al., 1977;Eckert et ai., 1973;Kamenski et al., 1987;Pirard and Kalitventzeff, 1978;Wong et al., 1974) in spite of the great advantages that they offer. The main reasons for this are the experimental difficulties and the requirements for high precision of the experimental data. The present paper is a critical analysis of the difficulties arising during the realization of this approach in laboratory practice. An approach is proposed that allows an effective application of the theory to the analysis of the experimental data. This approach has been worked out during the processing of the experimental data on ethylene oxide oxidation over supported silver catalyst. These data and
the derivation of the kinetic equation have been published recently (Petrov et al., 1988). Computer Programs The data processing was carried out by means of special programs in BASIC on a 9845B HewlettrPackard computer. These programs perform the following calculations: CHROM Program for Primary Kinetic Data Treatment. From the compositions, determined chromatographically, CHROM computes the following quantities: the partial pressures Pi (atm) of all the reactants and products both in the feedstock and in the converted gas mixture; the material balance of each experiment with respect to carbon, hydrogen, and oxygen; the CzHIO and 0,conversion degrees ( X Eand ~ Xo in %); the volume change for every experiment; the rates of consumption of the reactants (C2H40and 0,)and the rates of formation of the products (CO, and H20) W j (mol.h-'.(g of cat)-') by the formulas Wi = XiF"i/G (1) (2) Wco, = WH~O = ~XE#EO/G where X , is the degree of conversion of the corresponding reactant i (i = 1 is C2H40and i = 2 is O,), P i is the corresponding inlet flow rate (mol-h-'), and G is the catalyst amount in the reactor in grams.
0888-5885/91/2630-0639$02.50/00 1991 American Chemical Society
640 Ind. Eng. Chem. Res., Vol. 30, No. 4, 1991
A standard program "Basic statistic and data manipulation", supplied by the Hewlett-Packard Utility Library, was used for the statistical treatment of the experimental data. LINF Program. This program computes the rate R( 1) along the only independent reaction route-combustion of C2H40:
corresponding to a normal absolute error distribution and N
E2 = C(R(1)i" - R(l)i")2/(R(l)i")2
(8)
i=l
from the rates for reagents Wi, already determined by CHROM, using the relationships Wo, = -2.5R(1)
corresponding to a normal relative error distribution. Here N is the total number of experimental runs at all reaction temperatures, i is the number of the experimental run,the superscript m marks the rate estimated by the model, and the Superscript e denote the experimentally measured rate of C2H40 combustion. The NEM program calculatdd also the model deviation for each tested equation by the formula:
WE, = -R(1)
Md = (C@(l)im - R(l)rl/R(l)f)(100/M
(1)
CzH4O + 2.502 = 2C02 + 2H20
(3)
N
i-1
WH*O= 2R(1)
(4)
The system (4) consists of four equations with one unknown, R(1) ( W iare determined experimentally), thus being preset. Owing to the experimental error in the determination of Wi,this system becomes inconsistent. That is why we have to find its approximate solution. Overdetermined systems of linear equations can be solved by means of the Chebyshev approximation. The problem is equivalent mathematically to the linear programming problem to find min d within some preset limitations (Spivak et al., 1970; Andrushkevich et al., 1970; Slinko et al., 1972; Zhukhovitskii and Avdeeva, 1968): IPiR(1) - Wil I d
(5)
where Bi are the corresponding stoichiometric coefficients. It has been shown that a relative stability of the solutions exists at small experimental error of Wi (Spivak et al., 1974). The application of the standard programs for linear programming requires that the Haar condition is observed, which means that all the main minors of the stoichiometric matrix should necessarily be different from zero. This condition, however, is not always valid, and that is why we used the algorithm of Abdelmalek (Abdelmalek, 1975), which holds good also for the cases when the Haar condition is not observed. The experiments with greater error give nonfeasible solutions, and thus an additional check on the quality of the experiments is carried out. One of the authors of the present paper (Petrov and Shopov, 1977) has studied application of different computing procedures to determination of the rates along routes. It has been shown that in the cases when the law of error distribution is not known, the best result is achieved by using nonlinear programming methods with a minimization criterion:
E' = min max IWidc
- Wp'I/I
WpI
(6)
When the experimental errors are small, the linear and nonlinear programming give the same result. NEM Program. By means of the NEM program a kinetic equation fitting best to the experimental data is selected and the values of the preexponent and the activation energy of each kinetic constant as well as the reaction orders are estimated. The NEM program makes use of the partial pressures (determined by CHROM) and the rates along routes (determined by LINF). It is based on the NelderMead algorithm (Nelder and Mead, 1965). The following minimization criteria were used in the calculations: N
El = C(R(1)i" - R(l)i*)2 i31
(7)
(9)
The confidence intervals of the parameters were estimated by the Marquardt procedure (Marquardt, 1963). DISCRIM Program. We used the generalized criterion C (Hill et al., 1968) in order to determine simultaneously which model is the most probable and to estimate the values of its parameters:
C = WID + W2E
(10)
where the relative value of the criterion D is a measure of the divergence of the kinetic models and serves to discriminate between them. The relative value of the determinant E is a measure for the estimation of the values of the kinetic constants. Wl and W2are weight coefficienta by means of which the components of the criterion C acquire definite weight:
D
1 m-1
=-
2
m
i=1 j=i+l
I
pi(n)p.(n) 1
w1
a : - a?
(a2
+ uj2)(a2 + a?) [Rim-
= [m(l - P,)/(m - l)]"
w, = 1 - w,
(13) (14)
The calculations were carried out in the following sequence. An initial experimental grid of q (we accepted q = 10) experiments was constructed. The criterion used in selecting the initial grid was similar to the one described by Fiala and Fialova (1987) and relied mainly on the experience of the operator as well as on the potentialities of the apparatus. The grid covered the whole range of feasible values of the independent variables (outlet partial pressures, catalyst bed temperature). Then we accepted that Ej- = D , = 1and that all the models had one and the same a priori probability. An estimation of the parameters for all models was carried out on the basis of the initial set of q experiments. Further, the matrix M(q)and its corresponding determinant IM"J)Iand then D(q)were calculated and accepted to be E.,, and D,, correspondingly (Nathanson and Saidef, 1985): Sj = lMjl (15) fii
(16)
Ind. Eng. Chem. Res., Vol. 30, No. 4, 1991 641
For the remaining N - q experiments of the complete data set the determinant lMjl was estimated by using the relationship IMi'"+')l = IMj(")t(l
+ fiTDfiwi)
(20)
Then the values of C,D, and E were calculated once again. A next experiment was selected from the remainder of the complete set having a maximum value of the criterion C. It was added to the initial set to form a secondary set of q + 1 experiments. Once again the equations (11)-(20) were solved, this time for the secondary set. Then the a posteriori probability of the models was calculated, by using the equations m
Pj("+l) = Pj(n)P(Mj/~)[ ~P(")P(M~/cY)]-' (21) j=l
Finally the weight coefficients W , and W2were calculated and a new experiment number q + 2 was searched for among the remainder of the complete data set during the second trial. This iterative procedure continued until the posteriori probability of one of the models reached a value higher than 0.99. The kinetic constants evaluation was achieved by means of the nonlinear gradient method of Marquardt (Tabata and Ito, 1975) minimizing the s u m of the squared residuals El (7). As the experiments had been carried out at different temperatures, a reparametrization (Frank-Kamenetskii, 1947; Kittrell, 1970; Agarwal and Brisk, 1985) of the kinetic constants was imperative with the aim to (1) take into account the effect of temperature on the reaction rate and (2) equalize the scale of the preexponential factor It, and activation energy E , in the process of searching of the optimum parameters. Although the correlation between the parameters was decreased by this reparametrization, the magnitude of the confidence intervals of the real parameters did not change (Meites et al., 1987). This procedure, based on sequential experimental design, gave us the possibility to discriminate between the models with the simultaneous parameter estimation. Experimental Section Apparatus. The reaction rates were measured by a glass flow circulation system (Petrov et al., 1985). Differential pressure regulators were used to fix the same pressure (1.5 atm) for the C2H40,02,and Ar fed into the reactor. The flow rates were monitored by a Matheson multiple mass flow controller, Model 8249 the deviation did not exceed f1% of the set flow rate. The gas mixture flowed first through a six-way sampling valve and its composition was analyzed by GC (and corrected if necessary) before it was fed into the reactor. The fixed bed reactor with a preheater contained 25 g of catalyst and was placed in an oven. An electronic thermoregulator (PIT-3, USSR), maintained the desired oven temperature constant: the deviation did not exceed f l "C. The electromagnetic four-valve piston pump ensured an intense circulation of the reaction mixture at a rate of 800 L.h-'. The converted reaction mixture passed through a second six-way sampling
,
valve so that a sample was taken for analysis of the outlet mixture composition. The circulation cycle and the valve were placed in a thermostat at 110 "C to avoid water vapor condensation. Gas Chromatographic Analysis. The analyses were carried out by means of a Tzvet 110 gas chromatograph (USSR), equipped with a theimal conductivity detector and a 2-m-long Porapack Q column. Column temperature was 115 "C and the carrier gas was Ar (flow rate 30 d-min-'). Detector response was specially calibrated for each of the reagents. The volumes of the sampling loops were 0.25 mL for the inlet sample and 2 mL for the outlet sample. A microprocessor system ISOTCHROM-1 (Bulgaria) processed automatically the GC data. Experimental Conditions. The catalytic oxidation of ethylene oxide was studied under steady-state conditions at atmospheric pressure and over a wide temperature range (100-300 "C) (Petrov et al., 1988). The investigated feeds had the following composition: (1) 16.7% 02,33.3% C2H40,50.0% AI (1:2);(2) 25.0% 02,25.0% C2H40,50.0% Ar (1:l); (3) 33.3% 02,16.7% C2H40,50.0% Ar (2:l). The ethylene oxide contact time T was varied from 158 to 962 H-(g of cat).mol-'. Under these conditions the C2H40 conversion degree was in the range 0.5-25.5%. The system comes to a steady state within 30 min after a change of the experimental conditions. This was established by the statistical method of Wald-Wolfowitz (Himmelblau, 1973). The proceeding of the reaction in the kinetic region was checked by the Corrigan criterion (Corrigan, 1955). Reproducibility of the Kinetic Measurements and Stability of the Catalyst Activity. The gas chromatographic analysis was performed with 1% relative error. Every 2 weeks a standard test was performed and the catalyst activity was measured. The deviation never exceeded 8%. Catalyst. The catalyst, used in this kinetic study of C2H40oxidation, was synthesized in our laboratory. It has the following physicochemical characteristics: 20 wt % Ag, supported on a-A1203,promoted by Mg, Ca, K, and Na additives; Si02, added as a binding substance; spheroidal pellets, 6 mm in diameter; specific total surface area 0.14 m2.(g of cat)-', determined by the BET method (adsorption of Kr); Ag specific surface area 0.05 m2.(g of cat)-', determined by O2chemisorption at 200 "C; total pore volume 0.06 cm3.(g of cat)-'; prevailing pore radius 200 A. Results and Discussion Experimental Design. The application of the sequential experimental design methods to model discrimination and kinetic constant estimation encounters a number of experimental difficulties as it requires conducting experiments under strictly determined conditions: 1. For a gradientless reactor the independent variables are the temperature and the outlet concentrations (or partial pressures) of the reagents, the reaction rate being a function of them. However, in performing a kinetic experiment we are in a position to control only the inlet reactor parameters (space velocities, concentrations, partial pressures, etc.) but not the outlet ones. This is a specific feature of reactors with ideal mixing such as Berty's internal circulation reactor (Berty, 1979) and Temkin's external circulation reactor (Temkin et al., 1950; Temkin, 1979). Further, a sequential design of experiments is carried out with respect to the independent variables. It is practically impossible to conduct an experiment in such a way that a preliminary required composition of the outlet mixture is obtained. It would be possible only if we knew the kinetic equation in advance, but in fact this equation
642 Ind. Eng. Chem. Res., Vol. 30, No. 4, 1991
is the aim of the study. In this sense we were not able to control the independent variables during the experimental runs, but this has nothing to do with the experimental error. A way out of this situation can be found by one of the following methods: 1.1. Construction of an intermediate polynomial description, as it has been proposed by Petrov (Petrov, 1975). Such a description can be constructed by using the central composition rotatable plan (or any other experimental plan) in order to find a functional dependence between the inlet independent variables and the outlet concentrations and rates. Polynomials of the following type are constructed: Ci
f,(P"i,J'Oi,T)
ri = f 2 P i P i , T )
This description requires a considerable amount of experimental work and data processing without gaining essential information. Nevertheless it allows us to predict, with satisfactory accuracy, the values of the required quantities (Le., the conditions for the next experiment) and therefore to discriminate over a continuous experimental data field (i.e., to perform the next experiment under any conditions required by the sequential design). In spite of the fact that we are able to predict the values of Ci and ri, their experimental realization could be restricted by the limited apparatus potentialities or by some other circumstances (for instance explosive mixture formation). Therefore, in some regions the discrimination is carried out practically over a continuous experimental data field while in others over a discrete field. 1.2. Selecting the necessary conditions by experimental guesswork (Titova and Kokovin, 1982). This method, however, does not lead to satisfactory results, and besides it also requires a great number of experiments, yielding little information. 1.3. Carrying out a certain number of experiments and creating a discrete data field, over which the sequential design and the model discrimination are conducted (Bibin and Popov, 1973). The disadvantage of this approach is that not always the optimum experimental conditions for the discrimination between competitive models are selected beforehand. In this case the computer program selects these optimum conditions in the discrete data field after the experiments have been conducted. In return for this, however, the experimental work is reduced considerably. As our experience showed, the discrimination is effective upon uniform filling up of the experimental data field. In our case this means that the experimental runs were conducted at 250,275, and 300 "C and at three levels of the outlet partial pressures of O2 and C2H40for each temperature. This approach leads to a significant reduction of the number of experimental runs compared to the other two and guarantees reliable results. That is why we chose this approach in our work. 2. It is necessary to juxtapose models under "critical" experimental conditions (where the rates, predicted by them, differ considerably) so that the discrimination methods are effective. Such "critical" conditions, however, may also be created by the presence of a great experimental error and may lead to completely wrong conclusions. That is why it is extremely important to decrease the experimental error, as much as possible, by detecting and removing the experiments having a great error (i.e., the outliers). This is achieved either by statistical manipulation of the experimental results (Mah, 1982) or by using our LINF program prior to discrimination. It leads to a real solution only when the experimental error is slight,
while in the case of great experimental error it provides no solution (Spivak et al., 19T4). 3. The discrimination between the models is effective only if the values of their parameters do not differ too much from the real ones. Upon simultaneous discrimination and estimation of constants the initial experimental data set is used, on the bagis of which a preliminary tentative evaluation of the model parameters is made. The preliminary evaluation of the constants is rather inaccurate, because of the small number of these experiments, and it is quite possible that a wrong decision can be made by the program. In order to overcome these difficulties an estimation of the parameters (preexponents k" and activation energies E,) is made, based on the complete experimental data set (obtained as explained in paragraph 1.3.), using the algorithm of Nelder-Mead (NEM program). These estimates are then set in the DISCRIM program as starting values. Thereupon a much more stable work of the discrimination algorithm is achieved and the program operation time is diminished substantially. The approach, used by us for simultaneous discrimination and parameter estimation, consists of the following stages: I. A certain number of experiments are carried out, which cover uniformly the whole range of experimental conditions of practical interest. These conditions are the outlet reagent concentrations and the temperature. According to literature data published so far, the discrimination is conducted on the basis of experimental results, obtained at one and the same temperature. This fact affects considerably the reliability of the model selected as most probable. In the discrimination procedure that we worked out, the temperature is considered as independent variable. This problem has been discussed by Agarwal and Brisk (1985). This led to certain complications in the program; however, the obtained result is more reliable. 11. A statistical manipulation of the complete experimental data set is conducted with the aim to detect and remove the outliers. 111. The parameter values of each model are determined by the NEM program, based on the complete data set. These calculations are used for preliminary indirect discrimination between the tested models. IV. Parameter estimates, thus obtained, are used as starting values in the DISCRIM program, by means of which a simultaneous discrimination between several of the best models, selected during the previous stages of data treatment, is conducted. V. Additional experiments are carried out to check the reliability of the selected kinetic model. Kinetic Models and Reaction Mechanism. A great number of theoretically possible mechanisms of C2H40 oxidation over Ag catalysts has been considered. Different versions have been worked out, supposing that the C2H40, adsorbed on the catalyst surface, is oxidized either by the gas-phase 02,or by dissociatively or molecularly adsorbed oxygen, or by both. CHBCHOand HCHO as adsorbed species have been considered as possible intermediates of C2H40 oxidation. Various assumptions have been made with respect to the reversibility of the different elementary steps and the correlation between their rates. As a result of this analysis, a large number of kinetic equations have been derived and they have all been tested at the description of our experimental data.
Ind. Eng. Chem. Res., Vol. 30, No. 4, 1991 643 Table I. Values of Preexwnential Factors. Activation Enerdee. and Their 95% Confidence Intervals ~~
preexponents do 16.1 4.4 x 10-10 9.91 x 10" 3.37 x 10-14 7.03 X
lower limit 15.4 4.26 X 9.54 x 10" 3.29 x 1 0 4 4 6.86 X lo-"
dl d2
d3
d4
4
p:
~
upper limit 16.9 4.54 x 10-10 1.03 X lo-' 3.45 x 1 0 4 4 7.21 X lo-"
~
lower limit 9.42 18.80 12.64 30.93 22.56
activation energies, kcal/mol Ea upper limit 9.57 9.71 19.42 20.04 12.91 13.18 32.39 33.85 23.42 24.29
I
h
STEP NUMBER Figure 1. The a posteriori probability of the four simultaneously discriminated kinetic models as a function of step number.
The following selected best models, corresponding to different mechanisms, have been discriminated simultaneously by this procedure:
R(1)" =
dlPE#O
(1
+ d2Po + d3P0"~ + d4PEo + dSpE,,P04.5)2
(26) The kinetic constants of model 23 are denoted by ai and those of models 24,25, and 26 by bi, ci,and di, respectively, to avoid misunderstanding. The calculations, carried out in accordance with the scheme described above, showed that 16 trials were sufficient for effective discrimination between the competitive models; Le., only 26 experiments from the complete data set were used, at this stage, for the calculation of the criterion C (Hill et al., 1968). After 16 trials the probability P4for model 26 reached a value close to 1 (Figure l),which was the criterion for termination of the calculations: p1= 1.44 x 10-3; p2 = 8.54 x 10-5; P3 = 2.31 X lo4; P4 = 0.998 3290 So the best description of our experimental data was provided by model 26. It corresponds to the following mechanism- (27) (Petrov et al., 1988): step 1: Z+02'ZO2 2 step 2: 22 + 02 2 2 0 1/2 step 3: Z + C2H40 F? ZC2H40 1 step 4: ZC2H40+ ZO 2ZCH20 1 step 5: ZCH2O + Z 0 2 COP + H2O + 22 2 net: C2H4O + 2.502 = 2C02 + 2H20 (27) where Z is a vacant site on the catalyst surface. According
'
--
100 303 500 700 900 no0 CONTACT TIME r (hg-cu/fno/)
Figure 2. Juxtaposition of model and experimentally measured C2H40conversion degree as a function of contact time. (a,experimental conversion Xme;0, model conversion XEomat 250 O C and feed composition 25% 02,25% C2HI0, 50% Ar).
230 250 270 290 310 330
TEMPERATURE, O C
Figure 3. Comparison of model and experimentally measured C2H40conversion degree versw the temperature. (0,experimental at contact time 158 h a conversion XEoe;0,model conversion XEOm of catsmol-l (1)and 473 h-g of catemol-' (2) and feed composition 25% 02,25% C2H10, 50% Ar).
to this mechanism the adsorption of CzH4O and O2 (both dissociative and nondissociative) are fast steps a t equilibrium. A hypothetical intermediate of ethylene oxide oxidation is adsorbed formaldehyde. It is formed and oxidized further to C 0 2 and H 2 0 with the participation of both atomic and molecular adsorbed oxygen (ZO and Z02,respectively). Steps 4 and 5 are slow and irreversible; therefore their rates determine the total rate of the process. We used the method of Temkin (Temkin, 1964) for deducing model 26 from mechanism 27. The model is valid for the temperature range 250-300 "C (Petrov et al., 1988). The values of the preexponential factors and the activation energies of the kinetic constants and their 95% confidence intervals are listed in Table I. Model Adequacy Illustration. Figure 2 compares the model and the experimental dependence of C2H40conversion degree XEOon the contact time: XEO= fl(T) a t 250 "C and feed composition 25% 02,25% C2H40,50% Ar (1:l). The model describes quite well the increase in XEO with the increase in the contact time 7. The model and experimental C2H40conversion degree as a function of the catalyst bed temperature are juxtaposed on Figure 3, XEo = f2(t"). Both experimental series are for feed composition 1:1, but the contact time is
644 Ind. Eng. Chem. Res., Vol. 30, No. 4, 1991
Ci = outlet concentration of reagent i
D = relative value of the criterion of Box and Hill
FEED RATIO 02:C2ff40
Figure 4. Juxtaposition of model and experimentally measured CzHIO conversion degree as a function of feed composition. (0, experimental conversion Xme; 0,model conversion Xmmat 250 O C and contact time 315 h.g of cat-mol-' (1) and at 275 OC and 473 h.g of catsmol-' (2)).
different-158 and 473 h-g of cat-mol-'. As can be seen, the congruence between model and experimental values is very good. The model adequacy at the description of the effect of the feed 02:C2H40ratio on the C2H40 conversion degree, XEO= f3(Co:Cm) feed, is illustrated on Figure 4. The experimental series at 250 O C and contact time 315 hog of cat-mol-' (1)and at 275 "C and 7 = 473 h-g of catomol-' (2) are represented. The correspondence of model and experimentally observed C2H40oxidation rates as functions of the outlet partial pressures of the reagents is also very good (Petrov et al., 1988). The average deviation of the model is 17%.
Conclusions 1. The application of the methods of sequential experimental design shows that they do not lead to essential decrease in the experimental work. At the same time they require higher precision of the experiment. The computing time is considerably increased. In contrast to these disadvantages the reliability of the selected kinetic model and its corresponding mechanism of the process is notably higher. 2. The methods of sequential experimental design operate effectively and reliably only upon discrimination of a small number of models. The selection of these models is carried out by preliminary sifting of the pussible models and their corresponding mechanisms on the basis of statistical and chemical considerations. Upon discrimination of a greater number of possible models, the operation of the program becomes unstable and often this leads to accidental selection of a model. 3. The most effective approach for conducting experiments is the creation of a discrete field of the experimental data. 4. The realization of the experimental study by means of sequential design in a plug-flow reactor is quite simpler than that in an ideal-mixing reactor. This is due to the fact that in the former case the independent variables are the inlet parameters, while in the latter case they are the outlet parameters. 5. The selection of a model is more reliable when the temperature is also considered as an independent variable. 6. The parameter values, estimated in advance, should not differ too much from the most probable values so that the sequential experimental design procedure is effective. Nomenclature C = generalized criterion of Box, Hill, and Wichern
E = relative value of the determinant of Fisher's information matrix E' = minimization criterion fi = partial derivative of the reaction model rate with respect to kinetic constant ki P i = inlet flow rate of reagent i G = catalyst amount, g m = number of models M = Fisher's inverse matrix Mi = Fisher's inverse matrix, referring to model j Md = model deviation n, n + 1 = trial numbers N = total number of experimental runs POi = inlet partial pressures Pi(n) = a priori probability of model j Po = oxygen outlet partial pressure PEo= ethylene oxide outlet partial pressure P, = maximum probability value, reached by the best model P ( M j / a )= distribution of the density of probabilities a for the model j R(l)im,R(l)jm = rates of models i and j R(l), = experimental rate ri = rates for reagents, related to the whole system ss("+')/u? = normalized sum of squared deviations Si= determinant of Fisher's information matrix for model to = catalyst bed temperature, "C XEO= ethylene oxide conversion degree X i = degree of conversion of reagent i w = coefficient W1, W , = weight coefficients W i= rates for reagents, related to 1 g of catalyst Greek Symbols pi = stoichiometric coefficient of reagent i p = coefficient u* = variance of the experiment ai2 = variance of model i T = ethylene oxide contact time
Literature Cited Abdelmalek, N. Chebyshev Solution of Overdetermined Systems of Linear Equations. BIT 1975,15, 117. Aganval, A. K.; Brisk, M. L. Sequential Experimental Design for Precise Parameter Estimation. Znd. Eng. Chem. Process Des. Dev. 1985,24, 203-207, 207-210. Andrushkevich, M. M.; Buyanov, R. A.; Timoshenko, V. I.; Spivak, S. I. Kinetics of Oxidative Dehydrogenation of n-Butenes Over a Cr-Ni-Ca-phosphate Catalyst. Kinet. Catal. 1970, 11, 1419 (in Russian). Berty, J. M. Testing Commercial Catalysts in Recycle Reactors. Catal. Rev. Sci. Eng. 1979, 20 (l), 75-96. Bibin, V. N. Automated System for Calculation of Kinetic Models of Catalytic Reactions. In Chim. Kinetics and Catal.: Probl. Makrokinet. Metod. Vopr.; Kiperman, S. L., Ed.; Otdel, Inst. Khim. Fiziki: Chernogolovka, USSR, 1987. Bibin, V. N.; Popov, B. I. Application of the Methods of Mathematical Statistics on the Kinetic Study of Methanol Oxidation on Ferro-Molybdena Oxide Catalyst. Ref. Zh. Khim. 1974; Partsial'noe Okislenie Uglevodorodov. Metod. Math. Vopr. Kinet. 1973, 56-60 (in Russian). Box, G. E. P.; Lucas, H. L. Design of Experiments in Nonlinear Situation. Biometrica 1959, 46, 77. Box, G. E. P.; Hunter, W. G. Sequential Design of Experiments for Nonlinear Models. University of Wisconsin, Dep. Stat. Tech. Rep. 1963, No. 21. Box, G. E. P.; Hunter, W. G. Sequential Design of Experiments for Nonlinear Models. Proc. ZBM Sci. Comput. Symp. Stat. 1966, 113.
Box,G. E. P.; Hill, W. J. Discrimination Among Mechanistic Models. Technometrics 1967, 9, 57.
Corrigan, T. E. Interpretation of Kinetic Data-I. Fundam. 1955, April, 199.
Chem. Eng.
Ind. Eng. Chem. Res., Vol. 30,No. 4, 1991 645 Dumez, F.; Froment, G. F. Dehydrogenation of 1-Butene into Butadiene. Kinetics, Catalyst Coking, and Reaction Design. Znd. Eng. Chem. Process Des. Deu. 1976,15,291. Dumez, F.; Hosten, L. H.; Froment, G. F. The Use of Discrimination in the Kinetic Study of 1-Butene Dehydrogenation. Znd. Eng. Chem. Process Des. Dev. 1977, 16, 298. Eckert, E.; Hlavacek, V.; Marek, M. Catalytic Oxidation of CO on Cu/Al2O9. I. Reaction Rate Model Discrimination. Chem. Eng. Commun. 1973,1, 89-94,95-102. Fiala, J.; Fialova, J. Experimental Design. Silikaty (Prague) 1987, 64, 545. Frank-Kamenetskii, D. A. Diffwionand Heat Transfer in Chemical Kinetics; Academy of Sciences of USSR Press: Moscow-Leningrad, 1947,p 234.Froment, G. F.; Hosten, L. H. Catalytic Kinetics: Modeling. In Catalysis Science and Technoloay: Anderson, J. R., Boudart, M., Eds.; Springer Verlag: Berlin, FRG, 1981;Vol. 2, p 99. Hill, W. J.; Hunter, W. G.; Wichern, D. W. A Joint Design Criterion for the Dual Problem of Model Discrimination and Parameter Estimation. Technometrics 1968, 10, 145. Himmelblau, D. Process Analysis by Statistical Metods; Mir: Moscow, USSR, 1973;p 149. Hosten, L. H.; Froment, G. F. Proceedings of the 4th International Symposium on Chemical Reaction Engineering Heidelberg, 1976; Dechema: Frankfurt, FRG, 1976. Hsiang, T. C.; Reilly, P. M. A Practical Method for Discrimination Among Mechanistic Models. Can. J. Chem. Eng. 1971,49, 865. Hunter, W. G.; Reiner, A. M. Design for Discriminating Between Two Rival Models. Technometrics 1965, 7, 307. Hunter, W. G.; Hill, W. J.; Henson, T. L. Designing Experiments for Precise Estimation of All or Some of the Constants in a Mechanistic Model. Can. J. Chem. Eng. 1969, 47, 76. Kafarov, V. V.; Pisarenko, V. N.; Zhukova, T. B.; Merzhanova, R. F.; Mascheva, L. A.; Kuz’michev, A. I. A Sequential Criterion for Verifying Statistical Hypothesis. Dokl. Acad. Sci. USSR 1974, 218, 386. Kamenski, D.; Dimitrov, S.; Petrov, L. A Sequential Discrimination Procedure for Chemical Modeling-Catalytic Hydrogenation of Benzoic Acid. Sixth International Symposium on Heterogeneous Catalysis; Shopov, D., Andreev, A., Palazov, A., Petrov, L., Eds.; Publishing House of the Bulgarian - Academy of Sciences: Sofia, 1987;Part 1, p 99. Kittrell, J. R. Mathematical Modeling of Chemical Reactions. In Advances in Chemical Engineering: Drew. T. B.., Ed.:, Academic Press: London, UK, 1970; p 8. Mah, R. S. Detection of Gross Error in Process Data. AZChE J. 1982,28, 828. Marquardt, D. An Algorithm for Least-Squares Estimation of Nonlinear Parameters. J . SOC.Znd. Appl. Math. 1963, 12, 431. Meites, L.; Fanelli, N.; Papoff, P. The Dependence of the Variances of the Parameters in Nonlinear Regresion Analysis on the Number of the Data Points. Anal. Chim. Acta 1987,200, 387. Nathanaon, M. H.; Saidel, G. M. Multiple-objective Criteria for Optimal Experimental Design. Am. J. Physiol. 1985,248, R378. Nelder, J. A.; Mead, R. A Simplex Method for Function Minimization. Comput. J. 1965, 7, 308. I.
Petrov, L. A. Kinetic Model of Oxidative Dehydrogenation of Butene.1 on Multicomponent Oxide Catalyst. Zfeterogeneow Catalysis, Proceedings of the Third International Symposium, Varna; Shopov, D., Andreev, A., Palazov, A., Petrov, L., Eds.; Publishing House of the Bulgarian Academy of Sciences: Sofia, 1975;p 310. Petrov, L.; Shopov, D. Calculations of Rates for the Various Routes of Complex Reactions. React. Kinet. &tal. Lett. 1977, 7, 261. Petrov, L.;Eliyas, A,; Shopov, D. A Kinetic Model of Steady State Ethylene Epoxidation Over a Supported Silver Catalyst. Appl. CQtQl.1985, 18,87. Petrov, L.; Eliyas, A.; Maximov, Chr.; Shopov, D. Ethylene Oxide Oxidation Over a Supported Silver Catalyst. I. Kinetics of Uninhibited Oxidation. Appl. C Q t d 1988, 41, 23. Pirard, J. P.; Kalitventzeff, B. The Use of Experimental Design to Study the Kinetics of Hydrogenation of Ethylene on CopperMagnesia Catalyst. Znd. Eng. Chem. Fundam. 1978, 17,11. Reilly, P. M. Statistical Methods in Model Discrimination. Can. J. Chem. Eng. 1970,48, 168. Reilly, P. M.; Blau, G. E. The Use of Statistical Methods to Build Mathematical Models of Chemical Reacting Systems. Can. J . Chem. Eng. 1974,52, 289. Rozycki, C. Statistical Processing of Analytical Results; A Microcomputer program. Chem. Anal. (Warsaw) 1987, 32,569. Slinko, M. G.; Spivak, S. I.; Timoshenko, V. I. On the Criteria for Determination the Kinetic Model’s Parameters. Kinet. Catal. 1972, 13,1570 (in Russian). Spivak, S. I.; Timoshenko, V. I.; Slinko, M. G. Application of Chebishev Method of Equalization to Constructing the Kinetic Model of Complex Reaction. Dokl. Acad. Sci. USSR 1970, 192, 580. Spivak, S. I.; Slinko, M. G.; Timoshenko, V. I. Estimation of Significance of Measurement Influence on Kinetic Models of Chemical Reactions. React. Kinet. Catal. Lett. 1974, I, 99. Tabata, T.;Ito, R. Effective Treatment of the Interpolation Factor in Marquardt’s Nonlinear Least-Squares Fit Algorithm. Comput. J . 1975, 18, 250. Temkin, M. I. Scientific Foundations of Production and Selection of Catalysts; Boreskov, G. K., Ed.; Siberian Branch of the Academy of Sciences, USSR Novosibirsk, 1964. Temkin, M. I. Kinetics of Steady State Reactions. Adv. Catal. 1979, 28, 173. Temkin, M. I.; Kiperman, S. L.; Lukyanova, L. I. Flow-Circulation Method for Kinetic Study of Heterogeneous Catalytic Reactions. Dokl. Acad. Sci. USSR 1950, 74, 763. Titova, E. E.; Kokovin, G. A. Application of Nonstatistical Aproach to a Nonlinear Inverse Problem. Math. Methody Khim. Thermodyn. 1982, 170 (in Russian). Wong, C.; Scharer, J. M.; Reilly, P. M. A Discrimination Among Microbial Iron Oxidation Mechanisms. Can. J. Chem. Eng. 1974, 52, 645. Zhukhovitskii, S. I.; Avdeeva, L. I. Linear and Convex Programming; Nauka: Moscow, USSR, 1968;p 207.
Received for review January 2, 1990 Revised manuscript received September 13, 1990 Accepted October 23, 1990
