Ind. Eng. Chem. Res. 1996,34, 483-487
483
Modeling Transient Tracing in Plug-Flow Reactors: A Case Study Eric Walter Laboratoire des Signaux et Systemes, CNRS-ESE, 91192 Gif-sur-Yvette,France
Luc Pronzato Laboratoire Znformatique, Signaux et Systemes, URA CNRS 1376, Sophia-Antipolis, 06560 Valbonne, France
Yee Soong Department of Energy, Pittsburgh Energy Technology Center, P.O. Box 10940, Pittsburgh, Pennsylvania 15236-0940
Masood Otarod Department of Mathematics, University of Scranton, Scranton, Pennsylvania 18510
John Happel* Department of Chemical Engineering, Material Science and Mining, Columbia University, New York, New York 10027
For heterogeneous catalysis, ordinary differential equations in time can be used to interpret transient tracer data obtained in a continuous stirred tank reactor or gradientless reactor in which a steady-state reaction is occurring. The same procedure is often used for plug-flow reactor modeling. However, for transient tracing data obtained with a plug-flow reactor operating under differential conversion, it is necessary to use a set of partial differential equations involving distance through the reactor as well as time. In this paper, we present a case study of methanation on a nickel catalyst to illustrate the modeling procedure needed to estimate the desired parameters.
1. Introduction It is often possible to obtain important information about heterogeneous catalytic reactions by the use of transient tracing in which isotopically marked species are employed (Happel, 1986). The same two basic types of apparatus are used as is the case in studies without tracers, namely continuous stirred tank reactors (CSTR) and plug-flow reactors. In principle, it is simpler to correlate tracer data obtained from a CSTR, because they correspond to a set of linear ordinary differential equations. In contrast, correlating the data obtained from a plug-flow reactor involves solving partial differential equations. However, the apparent simplicity of the CSTR procedure is balanced by the difficulty of using this type of reactor, with the necessity of ensuring a proper gradientless operation. Consquently, many investigations using tracers have employed plug-flow reactors. For example, several recent studies of the methanation reaction have been conducted using the latter (Winslow and Bell, 1985; Soong et al., 1986; Zhang and Biloen, 1986; de Pontes et al., 1987; Stockwell et al., 1988; Yokomizo and Bell, 1989; Krishna and Bell, 1991, 1993). In these studies, instead of solving the appropriate partial differential equations, the authors assumed the validity of pool-type models similar to those used for CSTRs. We have shown how to obtain and solve partial differential equations required to describe tracer data obtained with plug-flow reactors (Happel, 1986; Happel et al., 1990; Otarod et al., 1992). In this paper, we wish to demonstrate the applicability of this methodology to a case study involving methanation. The paper is
* Author t o whom correspondence should be addressed. 0888-5885/95/2634-0483$09.00l0
organized as follows. Section 2 explains why we selected the data by Soong et al. (1986). Section 3 describes the model structure to be considered, the identifiability of which is studied in section 4. Section 5 recalls the experimental data to be modeled in Section 6. We hope that this case study may encourage other researchers to follow a methodology using the partial differential equations required for plug-flow systems.
2. Design of Plug-Flow Experiments There are a number of factors that must be considered in the design of plug-flow experiments. The speed of response of analytical equipment such as a mass spectrometer as well as its accuracy must be taken into account. Methods of introduction of traced species and removal of samples should avoid causing abrupt changes in system pressure. Dilution of the reacting species by an inert carrier gas may be desirable to lengthen response time. It is necessary for the ratio of the bed length t o particle diameter to be sufficiently large to minimize longitudinal diffusion. Furthermore, the velocity of gas through the catalyst bed must be slow enough for the reaction to take place t o an appreciable extent by diffusion of reactants into the catalyst particles before the gas leaves the reactor. The conventional concept of effectiveness factor (Satterfield, 1980) is not appropriate to quantify this effect, because it refers to catalyst performance under steady-state conditions. A suitable approximation is provided by Crank’s model for unsteady state diffusion in a sphere (Crank, 1957) in order to compute the ratio between the time it takes for a molecule to pass through the bed compared to the time to reach steady state conditions. This makes it impossible to use a pool model to interpret data from
0 1995 American Chemical Society
484 Ind. Eng. Chem. Res., Vol. 34, No. 2, 1995
expressed by a linear approximation as conversion occurs, so that the following relationships apply:
Figure 1. Structure of the model.
a plug-flow system. To satisfy these various constraints, some compromises in reaction temperature and pressure may be necessary, but usually it is possible t o operate much closer to normal operating conditions than is often required in using more sophisticated surface science techniques. Soong et a1.k data (1986) appear to be especially appropriate in satisfylng the above constraints. Moreover, their study involved aging of a methanation catalyst, as contrasted with other studies in which only freshly regenerated catalyst was employed, and thus supplied additional useful information on catalyst properties. For these reasons, we have chosen their data t o demonstrate the application of the procedure that we have developed to a case in which actual experimental data are modeled.
3. Development of the Model Structure The model structure to be considered is described by Figure 1. It refers solely to the reaction path followed by atomic carbon, which is the species serving as a tracer. Of course, to describe the entire reaction system, it would be necessary to consider the simultaneous oxygen and hydrogen transfers that occur during methanation. The reaction intermediates Ni-1 and Ni-2 correspond to two different types of catalytic sites which may evolve with time on stream. We do not imply knowledge of the exact nature of each of these sites. Other model structures were considered in (Walter et al., 1986) for a different nickel-based catalyst using a CSTR. The model structure used here is similar to one of those considered by Soong et al. and was found suitable for data correlation. Our primary purpose here is to derive a typical appropriate partial differential equation model, as contrasted with the ordinary differential equation models used in previous studies with CSTR models, and to estimate its parameters from the experimental data. Transposition to other plug-flow models would be straightforward. It is assumed that no longitudinal mixing occurs in the reactor. Reacting gases enter the catalyst bed at w = 0 and leave it at w = W. Following the notation given in the Nomenclature, material balance equations over a differential section can be written as
FC0(w) = F - VW, CCo(w)= FCo(w)/F
(5)
p H 4 ( w )= Vw, CCH4(w)= PH4(w)/F
(6)
For a unit step up from zero tracer initial conditions, the boundary conditions are as follows: zCo(o,w) = 0, z ~ ~ ~ ( o=, 0, w ) zNi-l(o,w)= 0, zNi-2(0,W)0, zCO(t,O)= u(t) (7) where u ( t ) is the Heaviside step function. For a unit step down in tracer in a system initially saturated with tracer, they become zCo(o,w) = 1, zCH4(0,w)= I, z”-~(o,w) = 1, zNiS2(0,~) = 1, zC0(t,O) = 1 - ~ ( t(8) ) For any given experimental run, corresponding to a given time on stream, CNi-land CNi-2and the contributions u 1 and u2 of the two branches to the overall methanation velocity V are assumed to be constants. Following the method given in Happel et al. (1990) and Otarod et al. (19921, one then easily obtains the solution for zc& at the exit of the reactor for a unit step down as
where z = @IF is the residence time, and u(t-z) is a delayed Heaviside function, equal to 0 if t < z and to 1 otherwise. Since V = u 1 + u2 is known from independent measurements, u1 and u2 need not be estimated individually. The parameters to be estimated may then be chosen as 8 =(cN~-~,cN~-~,~~)T
(10)
together with the residence time z, which is only known approximately. It seems worth noting that the same type of two-exponential response is obtained as in eq 25 in Soong et a1.k paper. An advantage of our approach is that the coefficients of the sum of exponentials are now obtained as functions of the model parameters, namely the surface concentrations of the intermediate species and the contribution u 1 of the upper branch in Figure 1 to the overall methanation velocity V. 4. Structural Properties of the Model
+
ulzNi-l(t,w) upNi-2(t,W)(4) For a plug-flow reactor in differential mode, it can be assumed that the change in concentration can be
When parameters of physical importance are t o be estimated on the basis of experimental data, it is desirable to check that results of data correlation have more than a purely empirical significance. The use of nonlinear least squares or any other statistical curve fitting procedure does not furnish insight into the fundamental validity of the parameters obtained. For this purpose, a mathematical study of the properties of the structure of the model is required. If it turns out that there are several possible values of the parameters that correspond to exactly the same behavior, it will then be impossible to discriminate between them on the
Ind. Eng. Chem. Res., Vol. 34, No. 2, 1995 485 basis of the tracer data alone. Such questions of uniqueness of the parameters associated with a given model behavior correspond to the notion of structural identifiability (see, e.g., Walter (1982, 1987)). The problem of isotopic assessment of fundamental catalytic mechanisms by kinetic modeling and its relation with the notion of identifiability were considered in detail in a previous paper (Happel et al., 1986). In that paper, we especially considered interpreting data from tests conducted using CSTRs. In the present paper, we extend this treatment to plug flow reactors. The study of the structural identifiability of the model described by FiguTe 1 is particularly simple. A model with parameters 8 will have exactly the same behavior on a unit step down as a model with parameters 8* if and only if both solutions correspond to exactly the same some of two exponentials, which can be written
zCH4(t,W = A exp(s,(t-t))
o'06fi0
50
150
200
250
300
350
400
Figure 2. Effect of aging on CO conversion (based on data in Soong et al. (1986)).
of CO conversion 100
(11)
and
100
Time on stream (hours)
V=(T
Identifying the coefficients, we get the two following solutions for 8 for each run (corresponding to a given time on stream)
-1
\
0.075.
o
+
(1- A ) exp(s,(t-z))u(t-t)
0.1
2011 0.13 mL/(min-g) (NTP) (14)
where the percentage of CO conversion for any given time on stream is obtained from Figure 2. Symbols on Figure 3 present the results of transient tracing experiments at various degrees of catalyst aging. The values given on the ordinates correspond to the fractional marking of methane at the exit of the reactor, zC&(t,W). 6. Results
For each run,a set of parameter values was estimated by minimizing the unweighted least squares criterion
(13) All parameters are therefore only identifiable to a limited extent, in that each of them can take two possible values. Such a situation, where the parameters can only take a finite number of values, is described as local identifiability. Further resolution may be obtained upon consideration of the entire series of runs. In what follows, we tentatively interpret Ni-1 as the fresh catalyst and associate it with the sites with the largest specific surface activity. That is, we impose ul/CNi-l > u ~ / C With ~ ~ -this ~ .choice, the model parameters become uniquely identifiable. Note that the same type of identifiability problem would have been encountered with a CSTR model with two paths in parallel. 5. Experimental Data The design of the experimental equipment employed and the data to be used in this investigation are taken from the original study by Soong et al. (1986). For the convenience of the reader, the information pertinent to the present study is summarized in this section. Figure 2 presents the effect of aging on the percentage of CO conversion. The rate of C 0 2 production was negligible (less than 5% of the methanation rate), so that it was not necessary to model it. The ratio of the length of the bed to the diameter of the catalyst particles was larger than 300, and therefore the effective longitudinal diffusion could be neglected in the plug-flow model that we employed. The volume of catalyst is approximately 0.26 mL, which corresponds to W = 0.13 g. Flow rates of 20 mumin (NTP) of CO and 40 mumin (NTP) of Hz were used (NTP means 20 "C and 1atm). The overall conversion of carbon monoxide to methane is therefore given by
(15) where y(ti) is the measurement obtained at time ti and ym(ti,8,t)is the corresponding model solution, given by (9). It was assumed that the catalyst properties remained constant for a single run. The parameters estimated were 8 and the residence time t (which is better evaluated from the complete curve than from the initial behavior). In each run, the estimated value for t was found equal to a few seconds. Figure 3 presents the model responses thus obtained, together with the original data. Note that the ordinates correspond to fractional marking and not to the degree of conversion, which is given in Figure 2. Each run is associated with a different degree of conversion. Figure 4 presents the and CNi-2, and Figure 5 gives estimated values for CNi-l the estimated values of u1 and UZ, together with V. On Figures 4 and 5, approximate 20-confidence intervals are indicated by vertical bars for each parameter estimate, under the simplifying assumption that V is known without any error. These were obtained from the computation of the inverse of the Fisher information matrix under the hypothesis of a measurement noise corresponding to independently identically distributed zero-mean Gaussian variables (see, e.g., Ljung (1987)). Note that with the hypothesis made in section 4 that u1/CNi-l > u ~ C " -no ~ ,crossing of parameters for Ni-1 and Ni-2 occurred, so that this hypothesis is consistent with a smooth variation of the parameters with time on stream. As Figure 2 shows, catalyst activity decreases with age but approaches a constant value after 200 h on stream. The major proportion of activity is associated with the species Ni-1 corresponding to u1 on Figure 5.
486 Ind. Eng. Chem. Res., Vol. 34, No. 2, 1995
Preexpasure: 2000 s
I
+:
i
I
Time on stream:
+
313 h
i
0.41
0.2
400
200
600
800
lD00
1200
1400
ldoo
1800
Time on stream (hours)
Time (seconds)
Figure 3. Comparison of data obtained by Soong et al. (1986) with responses calculated for a model with the structure described by Figure 1.
Figure 6. Specific surface activity (rate of production of methane as a per unit of surface concentration) of Ni-1 (+) and Ni-2 (0) function of time on stream. (Points are plotted on the vertical axis at zero time.)
roughly modeled, because the time on stream of 0.5 h for the first tracer run is close to the time required to conduct a tracing experiment, during which time the catalyst performance is assumed to remain constant according to eq 1. This problem is not encountered during the remaining tracer runs when catalyst activity change is much slower. 7. Discussion
i
100
50
150
200
250
300
350
Time on stream (hours)
Figure 4. Estimated values of CNi-l( x ) and CNi-2(0) for various times on stream. (Points are plotted on the vertical axis at zero time.) Vertical bars correspond to approximate +2a uncertainty intervals.
-a 3 .P
: s
". . 0.09- * 0.08-
0.07~ 0.06i
0.04-
2 0.03 !
j
1
0.05 -
0.02 -
~
4
9 t
0.01-
od
50
1W
1SO
250
200
300
1
1 1I
3SO
Time on stream (hours)
Figure 5. Estimated values of u1 ( x ) and u2 (0) for various times on stream and associated value of V = u1 u2 (*). (Points are plotted on the vertical axis at zero time.) Vertical bars correspond to approximate f 2 a uncertainty intervals.
+
Its specific surface activity (rate of production of methane per unit of surface concentration) is much larger than that of the Ni-2 component (Figure 6). It seems therefore reasonable to identify the Ni-1 component with the fresh catalyst. Table 1 summarizes the best values found for the parameters of interest for each run. Also indicated are the corresponding values of the criterion (15). The very small values obtained for this criterion confirm that the analytical solution (9) is indeed a very good representation of the data. The Ni-2 species seems to be rapidly formed during the initial 16 h on stream. This increase can only be
The primary object of the present study was to develop the appropriate partial differential equations to model the behavior of plug-flow catalytic reactor systems by means of tracers and to estimate the corresponding parameters from experimental data. In the model employed, it was assumed that the catalyst properties remain constant during the period of tracer transient observation, which has also been the case in a number of studies reported in the literature (Happel et al., 1980; Otarod et al., 1983; Winslow and Bell, 1985; de Pontes et al., 1987; Stockwell et al., 1988; Yokomizo and Bell, 1989; Krishna and Bell, 1991,1993). For this case study, we chose to model the data of Soong et al. which are of special interest in that they were obtained by specific experiments conducted separately over a substantial time period, during which catalyst aging occurred. We have shown that a model based on partial differential equations describing the behavior of plugflow reactor systems can furnish insight into the nature of the surface nickel species involved. These equations were simple enough to allow an analytical solution, which considerably speeded up the estimation of its parameters. More complicated models would be required if it were desired to elaborate the nature of intermediates, which might entail the use of numerical solvers. A recent study of Yadav and Rinker (1992) using step-response kinetics of methanation without tracer in a CSTR reactor confirms the importance of the COM2 feed ratio on the composition of such species. In Happel et al. (19801, a model was suggested to describe the methanation of synthesis gas in a CSTR reactor, which at low rates of C02 production reduces to the single path corresponding to the fresh catalyst in our plug-flow model (Figure 1). As already mentioned, the data used in our present case study (Soong et al., 1986) involved a negligible rate of CO2 production, so that the basic model structure was the same in both cases. Problems of identifiability and distinguishability (see, e.g., Walter (1982,1987))are of paramount importance
Ind. Eng. Chem. Res., Vol. 34, No. 2, 1995 487 Table 1. Estimated Parameters and Associated Value of the Least-SquaresCriterion for Each Run time on stream 0.5 16 40 64 313
V1
0.0682 0.0645 0.0569 0.0518 0.0488
v2 0.0316 0.0268 0.0264 0.0269 0.0202
V
CNi-1
CNi-2
0.0997 0.0914 0.0833 0.0787 0.0690
6.9366 8.7894 11.1672 9.6657 11.7286
17.1251 29.5471 32.0479 34.0833 35.1310
when attempting to estimate the parameters of phenomenological models such as those considered here. I t would be of considerable further interest to extend this research into developing a model that would simultaneously take into account catalyst performance at any one time and the effect of catalyst aging with time on pertinent parameters.
Acknowledgment We are grateful for support of this project by NATO Cooperative Research Grant No. CRG 921123.
Nomenclature Cco = concentration of CO in gas phase, volume fraction C c b = concentration of CHq in gas phase, volume fraction CNi-l, CNi-2= concentrations of hydrocarbon-containing intermediates on the sites Ni-1 and Ni-2 per unit mass of catalyst, mL(NTP)/g F o ( w )= flow rate of CO at longitudinal position w in bed, mL(NTP)/s F H 4 ( w )= flow rate of CHI at longitudinal position w in bed, mL(NTP)/s F = initial flow rate of CO at w = 0, mL(NTP)/s j ( 8 , t ) = unweighted least-squares criterion t = time, s (.IT = transposition operator (used to transform row vectors into column vectors) u(t-z) = unit step function (Heaviside step) with value of 0 at t < t and unity at t 2 t V = overall reaction velocity per unit mass of catalyst, mL(NTP)/s.g) ui = step velocity per unit mass of catalyst, i = 1, 2, mL(NTP)/(s*g) w = longitudinal position in bed, expressed as the mass of catalyst located between this position and the inlet, g W = total mass of catalyst, g y(ti) = datum collected at time ti ym(ti,8,r)= model output associated with y(ti) zi(t,w)= fractional tracer marking in circulating stream of product i as variant with time and space, i = CO, Ni-1, Ni-2, CHq ,8 = dead space including voids in catalyst bed, mL(NTP) 8 = column vector of the parameters to be estimated t = residence time BIF, s 0 = estimated standard deviation for a parameter estimate
Literature Cited Crank, J. Diffusion in a Sphere. The Mathematics of Diffusion; Clarendon Press: Oxford, UK, 1957;Chapter 6. Happel, J. Isotopic Assessment of Heterogeneous Catalysis; Academic: Orlando, FL, 1986. Happel, J.; Suzuki, I.; Kokayeff, P.; Fthenakis, V. Multiple Isotope Tracing of Methanation over Nickel Catalyst. J. Catal. 1980, 65, 59-77. Happel, J.; Walter, E.; Lecourtier, Y. Isotopic Assessment of Fundamental Catalytic Mechanisms by Kinetic Modeling. I n d . Eng. Chem. Fundam. 1986,25,704-712. Happel, J.; Walter, E.; Lecourtier, Y. Modeling Transient Tracer Studies in Plug-Flow Reactors. J. Catal. 1990,123,12-20.
criterion 5.6 10-5 4.5 10-5 3.0 x 10-5 3.3 x 10-5 6.1 10-5
Krishna, K. R.; Bell, A. T. A n Isotopic Tracer Study of the Deactivation of RwTiOz Catalysts during Fisher-Tropsch Synthesis. J. Catal. 1991,130,597-610. Krishna, K. R.; Bell, A. T. Estimates of the Rate Coefficients for Chain Initiation, Propagation, an$.Termination During FischerTropsch Synthesis over R u t h e n i u m t a n i a . J. Catal. 1993, 139,104-118. Ljung, L. Parameter Estimation Methods. System Identification, Theory for the User;Prentice-Hall: Englewood Cliffs, NJ, 1987; Chapter 7. Otarod, M.; Ozawa, S.;Yin, F.; Chew, M.; Cheh, H. Y.; Happel, J. Multiple Isotope Tracing of Methanation over Nickel Catalyst, 111. Completion of 13C and D Tracing. J.Catal. 1983,84,156169. Otarod, M.; Lecourtier, Y.; Walter, E.; Happel, J. Mathematical Modeling of Transient Tracing in Plug-Flow Reactors. Chem. Eng. Commun. 1992,116,127-151. de Pontes, M.; Yokomizo, G. H.; Bell, A. T. A Novel Method for Analysing Transient Response Data Obtained in Isotopic Tracer Studies of Carbon Monoxide Hydrogenation. J. Catal. 1987, 104,147-155. Satterfield, C. N. Experimental Methods. Heterogeneous Catalysis in Practice; McGraw-Hill New York, 1980;Chapter 11. Soong, Y.; Krishna, K.; Biloen, P. Catalyst Aging Studied with Isotopic Transients: Methanation over Raney Nickel. J.Catal. 1986,97,330-343. Stockwell, D. M.;Chung, J. S.; Bennett, C. 0. A Transient Infrared and Isotopic Study of Methanation over NickeYAlumina. J. Catal. 1988,112,135-144. Walter, E. Identifiability of State Space Models; Springer-Verlag: Heidelberg, 1982. Walter, E., Ed. Identifiability of Parameteric Models; Pergamon: Oxford, 1987. Walter, E.; Lecourtier, Y.; Happel, J.; Kao Juo-Yu Identifiability and Distinguishability of Fundamental Parameters in Catalytic Methanation. AZChE J. 1986,32,1360-1366. Winslow, P.; Bell, A. T. Studies of the Surface Coverage of Unsupported Ruthenium by Carbon-and-Hydrogen-Containing Adspecies During Carbon Monoxide Hydrogenation. J. Catal. 1985,91,142-154. Yadav, R.; Rinker, R. G. Step-Responses Kinetics of Methanation over Nickel-Alumina Catalyst. Ind. Eng. Chem. Res. 1992,31, 502-508. Yokomizo, G. H.; Bell, A. T. Isotopic Tracer and NMR Studies of Carbonaceous Species Present During Carbon Monoxide Hydrogenation over R u t h e n i u d t a n i a . J.Catal. 1989,119,467482. Zhang, X.; Biloen, P. A Transient Kinetic Observation of Chain Growth in the Fischer-Tropsch Synthesis. J. Catal. 1986,98, 468-476. Received for review March 18,1994 Revised manuscript received September 28, 1994 Accepted October 5,1994@
IE9401723 Abstract published in Advance ACS Abstracts, December 15,1994. @
