5990
Ind. Eng. Chem. Res. 2002, 41, 5990-6004
Kinetics Modeling of Cellulose Fast Pyrolysis in a Flow Reactor D. Rizzardi Soravia and P. Canu* Dipartimento di Principi e Impianti di Ingegneria Chimica “I. Sorgato” (DIPIC), Universita` di Padova, Via Marzolo 9, 35131 Padova, Italy
We developed a model of nonisothermal cocurrent gas-solid flow reactor with the purpose of interpreting experimental data of cellulose fast pyrolysis obtained in our laboratories. Data were obtained in an entrained-flow reactor operated with a low solid-gas ratio. Avicell cellulose powder was carried by nitrogen into a quartz tube with different heating sections. Weight loss measurements were collected at temperatures between 400 and 650 °C, residence times between 50 and 140 ms, and nominal heating rates ranging from 6000 and 13 000 K/s. Proper modeling of these data includes accounting for the particle velocity and temperature profiles along the reactor, which in turn requires a detailed description of heat transfer inside the tube to predict the temperatures of the gas and radiant surfaces. The contribution of moisture release to the total weight loss has been included and quantified. The single-particle model without intraparticle temperature gradients has been used, because of the small particle size and solid flow rate compared with the carrier gas flow rate. Several reaction mechanisms are discussed, eventually demonstrating that only the inclusion of an activation step can explain all of the data, from the lowest to the highest temperature tested. Finally, we present a discussion of the actual existence of two regimes in cellulose pyrolysis, depending on the heating rate and temperature achieved, suggesting that this view could result from the limitations of using excessively simplified kinetics to interpret the experimental data. Introduction Fast pyrolysis of biomass is a useful process for converting a poor solid fuel into more valuable oil or gas fuels.1 The specification of “fast” is motivated by the high heating rates characteristic of the process. More insight can be obtained by splitting the heating rate concept into residence time and temperature achieved, where fast pyrolysis is always characterized by a short contact time and the temperature is not precisely specified. Accordingly, we observe different products depending on the temperature, with gas prevailing at higher temperature and oils abounding at low temperature. A short contact time remains a common requirement for either gas or oil to be obtained. It is still an open question whether the different distribution of products is a result of either physical processes or some sort of chemical mechanism or both. It is difficult to provide a definitive answer mostly because of the difficulties involved in the experimental investigation of a short-contact-time-high-temperature process involving a solid. We obtained some experimental data in an entrained-flow reactor (EFR) under carefully controlled conditions2 with the purpose of developing some advancement in the understanding of this process through modeling. We investigated the role of some frequently overlooked aspects, such as moisture removal, particle size, and aggregative flow, in addition to the more traditional, and relevant, effect of temperature. Here, we develop a model for such a nonisothermal reactor where the relevance of the physical transport process is discussed but the role of the chemical mechanism is eventually the most relevant information in
terms of explaining the experimental measurements. It is shown that a multistep reaction mechanism (or a semiglobal model, according to the classification of Di Blasi3) is a requirement for explaining the existence of two regimes in the process,4 which are apparently determined by the heating rate. We compared one onestep and several multistep mechanisms with the experimental data. Agreement with our experimental data and literature data obtained under comparable conditions5 calls for at least two steps, with the independent activation of cellulose preliminary to the degradation steps. One-step models can reasonably follow the experimental observations, but not all of them can do so with a single set of parameters. Many attempts have been made to attain a better description of the process using more complex schemes,6,7 involving distributed activation energy (DAE) models8 or deactivation models.9 Many authors10,11 agree on the ideas originally formulated by Broido et al.7 that the cellulose degradation begins with a mostly intermolecular process, a reduction of the degree of polymerization, with limited or null volatile production. The product of this initial stage, called “active cellulose” or “anhydrocellulose”, further reacts along two parallel (and thus competing) paths; one involves ring scission and gas and char as products, and the other is a depolymerization of the terminal groups in the chain, which produces tars. An example is reported by Bradbury et al.6 and is usually called the “Shafizadeh model” k1
cellulose 98 active cellulose k2
active cellulose 98 tar k3
active cellulose 98 0.35char + 0.65gas * Corresponding author. Tel.: (+39) 049-8275463. Fax: (+39) 049-8275461. E-mail:
[email protected].
The first step is very fast, especially at high heating
10.1021/ie020299p CCC: $22.00 © 2002 American Chemical Society Published on Web 10/23/2002
Ind. Eng. Chem. Res., Vol. 41, No. 24, 2002 5991
rates, and gives rise to little weight loss, so it has been often neglected, and its existence has even been questioned. The very low residence times provided by the EFR reactor used in this research allowed the importance of this step to be recognized, in particular as a cause for the two regimes observed in the literature.4 Many excellent and detailed reviews of cellulose pyrolysis are available in the literature;4,12,13 here, we discuss semiglobal models with reference to specific experimental data2,5 obtained at very short residence times over a wide range of temperature. The possibility offered by short residence times is unique in terms of the identification and characterization of fast steps in a reaction mechanism. Moreover, we try to demonstrate the need for the proper identification of the solid temperature in the discussion of any reaction mechanism. Because of the complexity of experiments at elevated heating rates, the prediction of the sample temperature is itself a significant contribution to the modeling. A similar study has, to our knowledge, only been attempted by Brown and co-workers,14 who developed detailed computations of particle trajectories, through computational fluid dynamics (CFD), and temperature, through a steady-state heat balance. Their scope was, in a sense, opposite to that addressed in this paper: they characterized their reactor and the particle behavior therein to confirm that their experimental results15 can be discussed under the assumption of a pure chemical regime, because they used small particles with diameters of 50 µm or less. Here, we use modeling of the particle temperature throughout the kinetic study, because we believe that, under fast pyrolysis conditions, the temperature history of the particles is extremely influential on the results but also depends on the reaction progress. Eventually, particle temperature and kinetics are always intimately coupled so that the two aspects must be studied together, and a distinction between the processes that determine the particle temperature and the kinetics is misleading. Di Blasi3 has thoroughly reviewed many semiglobal mechanisms that have been proposed in the literature and compared them in terms of weight loss kinetics and species selectivity, selecting different sets of parameters from the literature, but she made no reference to experimental data. The most likely multistep mechanism is also used in this work in a parametric study to demonstrate its consistency with both “slow” and “fast” pyrolysis experiments when the first step has a high activation energy. It is suggested that the dispute about the activation energy of cellulose pyrolysis might be a consequence of an excessively simplified kinetic model that reduces any chemical process to a single first-order reaction. Such a model can provide a good fitting of weight loss data but yields kinetic parameters that are useful only for a specific range of experimental conditions. This is due to the excessive simplification of the approach that nevertheless makes its success so extensive. Overview of the Devolatilization Model In the following discussion, we develop a model to describe and interpret the devolatilization process as carried out in our experiments in an entrained-flow reactor, the details of which are reported in the original paper.2 Several assumptions and specific developments are connected with the specific apparatus and conditions used. Loosely speaking, in the standard experiment, the
cellulose particles are carried into the reactor by an inert gas, and they are heated through convection with the carrier flowing along the hot reactor walls and radiation from the alumina tube. Particles can exchange momentum and heat among themselves through collisions. We neglected these contributions, assuming single isolated particles because the solid concentration in the experiments was always kept to the minimum compatible with acceptable experimental errors. Accordingly, the singleparticle model is representative of the entire swarm of particles, and interaction effects are neglected. Because most of the experimental results obtained in our work,2 and in the literature as well, are given in term of weight loss, the most relevant equations account for the variations of different components’ weights, that is
dWi
NR
)
dt
∑j νijRj - Ji
i ) 1, ..., NC
(1)
The material balances above include variations caused by the NR reactions, as well as diffusive fluxes from the solid to the gas phase, Ji. A stoichiometric matrix vij is required to connect the change for each species to any reaction that involves it. The sum of all species weights at any time t compared with the initial value provides the total weight loss during the process. Equations 1 are formulated in terms of time so that they can be compared to typical expressions reported in the literature to describe weight loss data in batch experiments [e.g., thermal gravimetric analysis (TGA)]. However, such a formulation is not suitable for describing EFR data, which are time-invariant during a single experiment. In a continuous-flow reactor processing a solid suspension, the relevant time is connected with a single particle, which is a function of its position in the reactor through the local velocity
dz ) vp(z) dt
(2)
so that the proper formulation of the material balances turns out to be
vp(z)
dWi dz
NR
)
∑j νijRj - Ji
i ) 1, ..., NC
(3)
It is immediately evident now that knowledge of the particle velocity along the reactor, vp(z), is crucial. At this point, it suffices to say that the particle velocity is tightly connected with the gas velocity, that changes with local temperature according to gas density variation, although the precise determination will be discussed later in this work. A second crucial issue not immediately clear from eqs 3 is that the reaction rates Rj vary dramatically along the reactor because of the variation of the particle temperature. In principle, particle temperature is a function of the position both along the reactor and inside the particle itself. As discussed below, small particles can be assumed to be isothermal, so that the solid temperature can be assumed to change with the position in the reactor only. Accordingly, the thermal energy balance of a single particle can be written as
dpFpcpvp dTp ) Q˙ conv + Q˙ rad + Q˙ evap + Q˙ R 6 dz
(4)
5992
Ind. Eng. Chem. Res., Vol. 41, No. 24, 2002
which includes several mechanisms that modify the temperature, namely, convection and radiation, phase change inside the particle, and heat of reaction. Each of these terms will be specified in the following discussion, when appropriate. Note that the heat of reaction depends on the reaction rates and enthalpies, thus coupling eqs 3 and 4, which must be solved together. The mass and energy balances in eqs 3 and 4 clearly state that the EFR is not isothermal, and a proper investigation on the kinetics must account for variation of the local temperature of the solid. Another important feature of the EFR is the residence time, which depends on the particle velocity, vp(z), again a local function that is related to the temperature profile of the gas. It can be calculated by the integration of eq 2 from the entrance of the reactor to the reactor length, Ztot
ϑ)
∫0Z
tot
dz vp(z)
mp
Gas and Particle Flow in the Reactor Significant simplifications are introduced by the plugflow and single-particle (no-aggregation) approximations. A brief discussion of these assumptions is mandatory. Plug flow for both the carrier gas and the particles has been precisely determined through local measurements in riser and downer tubes.16 These measurements demonstrated that downflow allows for the rapid stabilization of a plug-flow profile for the gas and the particles, unlike the upflow (riser) regime. The conditions in our reactor are somewhat different from the conditions in those experiments,16 but the differences allow us to be even more confident about the plugflow assumption. The inlet velocity in our reactor is the same as that used in the earlier study,16 but the velocity then increases to up to 3 times the initial value after the particles are heated, thus entering further in the turbulent regime where plug-like profiles are approached. The particle residence time distribution (RTD) was also measured,16 and it was demonstrated that a downer reactor produces a very narrow distribution after quite long distances, which is definitely not the case for upward flow. The single-particle assumption can be discussed by evaluating the average solid density in the reactor at the highest solid flow rate (1 g/min, which is indeed very seldom used) and the lowest velocity (4.24 m/s). Under these circumstances, we calculate a density of 8.4 × 108 particles/m3, corresponding to an average distance between particles of 1.06 × 10-3m, or 12 particle diameters. With the lowest solid feed rate used ( 100
() µb µw
0.14
2100 < Re < 104 Re > 104
(12)
5994
Ind. Eng. Chem. Res., Vol. 41, No. 24, 2002
Figure 4. Experimental (circles) and calculated (solid line) gas temperature profiles.
Figure 6. Temperature profiles (a) of the gas and of nonreactive cellulose particles in the cases of heating (b) by both convection and radiation and (c) by convection only. T set ) 500 and 900 °C.
alumina wall through the air and the quartz tube. The simplified form of the energy balance (eq 4) can be written as
Fpcpvp
Figure 5. Radial temperature profile (T set ) 500 °C at the center of the reaction tube, z ) 0.2 m).
taking into account the inlet effects through the Graetz number. TR1, the temperature of the inner side of the reactor walls, is a function of TR2, which depends on TR3, and so on until TR4. Accordingly, eq 11 must be solved in conjunction with the conservation of heat flux across each layer (three additional algebraic equations). The results are shown in Figure 4 for 400, 450, and 500 °C. The agreement is quite good, except for the cooling at the end of the reactor, which is not modeled. Calculation confirms that most of the temperature drop is concentrated in the interstice of air, as shown in Figure 5, so that the temperature of the inner walls is not much higher than the gas temperature, even in the middle of the most active heating element (the first). It is worth mentioning that the temperature reported by the thermocouple used to control the power supplied to the corresponding heating element that is placed in the air interstice agrees with the model predictions. Its reading of 590 °C is consistent with the temperature excursion of 565-610 °C in this region, given that the thermocouple is almost as thick as the space between the quartz and alumina. Now confident with the thermal model of the gas, we can use the predicted wall temperature as the basis for the evaluation of the radiant contribution to the particle temperature. An initial energy balance on the particles can be formulated, with the aim of determining the dynamics of the particle temperature along the reactor in the absence of any reaction or phase change with significant enthalpy change. A complete energy balance will be developed and used later, after the discussion of the reaction mechanism and other devolatilization processes. So far, we include two mechanisms of heat exchange for the flowing cellulose particles, namely, convection with the carrier gas and radiation from the
dTp Sp ) (Q˙ conv + Q˙ rad) ) dz Vp 6 [h (T - Tp) + σ(TR34 - Tp4)] (13) dp gp g
where the specific heat of cellulose has been estimated19 as ∼1.67 kJ/(kg K) and the emissivity of the solid20 as 0.95. For the convection coefficient, we used the standard correlation21
Nugp ) 2.0 + 0.60Rep1/2Pr1/3
(14)
Gas properties were calculated at the film temperature, i.e., at a mean temperature between Tg and Tp. The results are shown in Figure 6. We can immediately observe that the isothermal region of the reactor is even shorter for the particle, underlining the importance of proper modeling of the particle temperature history in the reactor in the discussion of any sort of kinetic model. Moreover, we can conclude that the role of radiant heat flux is rather limited at 500 °C, one of the highest temperature used in the experiments. Calculations extrapolated to 900 °C, also shown in Figure 6, indicate that, at such high temperatures, radiation is a key mechanism in rapidly heating the particle, avoiding the initial delay of the particle temperature typical of convection. Nevertheless, in the following calculations, radiation has always been kept in the model. The calculations of particle temperature obtained with eq 13, and shown in Figure 6, are fairly simplified and preliminarily reported here to support the need for a proper accounting of variations in Tp along the reactor to begin any kinetics discussion. More accurate formulations of the particle energy balance will be developed later, in conjunction with a discussion of the kinetics of cellulose degradation and moisture release. A few general comments on such an energy balance are in order. In the particle energy balance that we use, i.e., eq 13 and the more detailed formulations to come, we assumed a unique temperature for the entire particle, neglecting any internal temperature distribution. This point is rather critical and will be thoroughly discussed in a specific section in the following, when the quantitative results of the simulations are available.
Ind. Eng. Chem. Res., Vol. 41, No. 24, 2002 5995
Another critical approximation could be the constant value for the solid heat capacity, because, during devolatilization, the porosity of the solid is expected to change, thus affecting its properties (specific heat, but also density and thermal conductivity). Because we used small particles and rather limited conversions, we neglected this secondary effect, the quantitative modeling of which is severely hindered by the uncertainties in the physical properties and structure of the degrading particle. We also neglected any effect on the gas temperature due to the occurrence of the pyrolysis because of the very low flow rate of the solid compared to the gas, regardless of the reaction heat. Note that the profiles of the gas temperature, wall temperature, and velocities of both the gas and the particles were determined just once, with different solutions being obtained for the each nominal experimental temperature, T set. These results are not significantly affected by the pyrolysis behavior, but only by the temperature of the experiment. Accordingly, vp(z), vg(z), Tg(z), and TR1(z) were stored in the form of piecewise polynomials and used for the following kinetic study.
Figure 7. Weight loss due to drying (experimental and calculated).
The value of kH2O was determined with the results of the drying test carried out at 110 °C, where no reaction occurs. Equation 18 must be solved together with the energy balance (eq 4), which, in this case can be written as
dpFpcpvp dTp ) Q˙ conv + Q˙ rad + Q˙ evap 6 dz
Effect of Initial Moisture The water content in the cellulose powder can be significant, so that part of the weight loss is due to drying. The phenomenon must be isolated to avoid confusing it with a contribution from the reaction. Because specific experiments have been performed to ascertain the drying kinetics,2 we can develop a model based on those data. The initial moisture in our cellulose was about 5 wt %, the concentration at which the release by evaporation is controlled by internal diffusion of the liquid.18 The process can be described as firstorder in term of the solid moisture, χH2O ) WH2O/W 0DC
dχH2O dt
) -k′H2O(χH2O - χe)
(15)
where k′H2O is proportional to the internal diffusivity of liquid divided by the square of the particle diameter. Assuming that the equilibrium moisture content of the solid, χe, is zero (because the steam content of the carrier is null at the inlet and cannot increase significantly with this process because of the small solid concentration and because the diffusion coefficient is proportional to Tp), eq 15 can be reformulated as
dχH2O dt
) -(kH2OTp)χH2O
L H 2O ≡
- W H 2O
W 0DC
0 ) χH - χ H 2O 2O
(16)
(17)
so that eq 16, together with the result of eq 2, can be reformulated as follows
vp
dLH2O dz
0 ) -kH2OTp(χH - LH2O) 2O
with a nonnegligible contribution due to the evaporation process, accounting for the heat of vaporization of water and the rate of steam release. The additional term can be expressed through eq 18, resulting in
cp
dTp 6 [h (T - Tp) + σ(Tem4 - Tp4)] ) dz dpFpvp gp g dLH2O evap (T ) ∆HH (20) O p 2 dz
evap ∆H H was calculated as a function of particle tem2O perature at 1 bar, assuming a value of 2258 kJ/kg at 100 °C. Because the experimental data are given in terms of global weight loss, L, a connection to the moisture content χ must be formulated. The global weight loss is defined in terms of the actual and initial masses of solid material
L≡
W0-W W0
(21)
Specifying the contribution to the material mass of both water and dry cellulose, we obtain
For sake of consistency with the following equations, we introduce the weight of moisture loss, LH2O, defined as 0 WH 2O
(19)
(18)
WH2O + W 0DC W ) L)1- 0)1- 0 W W H2O + W 0DC 1-
WH2O + W 0DC 0 (χH + 1)W 0DC 2O
)
LH2O 0 1 + χH 2O
(22)
The results of data fitting are shown in Figure 7. The agreement is extremely good, and the parameter values estimated are
kH2O ) 0.040 1/(K s)
and
0 χH ) 0.048 2O
with the initial moisture content coincident with the value determined independently by oven dehydration.2
5996
Ind. Eng. Chem. Res., Vol. 41, No. 24, 2002
All of the following calculations of weight loss due to pyrolysis include this submodel of moisture removal, because it is very rapid, taking place as the first step of any devolatilization, and it contributes to the total weight loss with a nonnegligible value when considering low-conversion experiments. The removal of water always completes at values of temperature slightly above 100 °C, as suggested by Figure 7, regardless of any subsequent temperature rise. This removes any concern about limitations of the water diffusion model, with the diffusion coefficient linearly depending on T. Pyrolysis Kinetics: Pseudo-First-Order Kinetics A single first-order reaction (SFOR) is the most frequently used kinetic model for describing the pyrolysis of complex solid fuels22 and cellulose as well.12 A simple pseudo-first-order model assumes that the rate of solid consumption is proportional to the amount of remaining pyrolysable cellulose
dWDC ) -k(WDC - W∞DC) dt
(23)
The presence of W ∞DC is the main disadvantage of this model, because the final weight (unpyrolyzed residue) is not a property of the cellulose but rather depends on the “thermal history” of the sample. The parameters k and W ∞DC obtained by fitting eq 23 to experimental data are not expected to be useful in scaling-up the process, because they are strongly dependent on the kind of experiment and operating conditions used to obtain the data. We nevertheless decided to apply this model preliminarly in trying to explain our data. Equation 23 must be solved together with eq 18 for the rate of water loss, with eq 20 being solved for the temperature of the reacting particle. All of the differential equations must be formulated in terms of z, rather than t, through eq 2, calling for the particle velocity submodel results. Once the profiles WDC(z) and LΗ2Ο(z) are known, the total weight loss L to be compared with the experimental data can be calculated as
L)1-
WH2O + WDC LDC + LH2O W ) 1 ) (24) 0 W0 W0H2O + W 0DC 1 + χH 2O
using the weight loss relative to dry cellulose, defined as
LDC ≡ 1 -
WDC W 0DC
(25)
The integration of differential eqs 18, 20, and 23 was carried out numerically, using a variable-step-variableorder solver based on the numerical differentiation formulas, as implemented in MatLab’s ode15s function. Such an algorithm is suitable for stiff ordinary differential equations, which is frequently verified in multistep mechanisms such as those discussed in the following. For consistency, we always used the same numerical technique for all of the mechanisms. Because the resulting profiles of Tp, WDC, and χΗ2Ο depend on the values of the kinetic constants, we identified the optimum values by fitting L(z) from eq 24 to the experimental data at the available positions (zi ) 0.4, 0.8, 1.2). Because the standard deviation σi of the
Figure 8. Calculated and experimental weight loss at three temperature levels. Error bars indicate (σ intervals. Single firstorder kinetics with neglect of the reaction heat. The best-fit parameters are ln(A) ) 19.02 (A in s-1), E/R ) 12 860 K, L ∞dc ) 0.55.
experimental measurements was quantified through several replicas of the same experiments, we mimimized the weighted sum of error, defined as
∆)
x∑(
)
Lcalc,i - Lsper,i
i
σi
2
(26)
The results of fitting an SFOR model to our hightemperature data are shown in Figure 8. The estimated activation energy is 107 kJ/mol, smaller than the mean literature value4 for fast pyrolysis (140 kJ/mol), which, in turn, is much smaller than values obtained with fast TGA.12,13 The SFOR model could not fit the results obtained at 400 and 450 °C, because of their particular behavior, where we observe quite a rapid increase of conversion at the beginning of the reactor, followed by a very slow reaction progress. Note that the calculation reported in Figure 8 used a simplified form of the energy balance, namely, eq 20, where a zero heat of reaction was assumed. The question of whether cellulose pyrolysis is endothermic, exothermic, or both has frequently been discussed. The experimental indications are confusing because many processes are taking place during pyrolysis, some endothermic, others exothermic. Once again, some information required to identify the kinetics (i.e., the particle temperature, which depends on the reaction heat) depends on the actual mechanism, creating a short circuiting that does not allow the two problems to be separated. Moreover, it turns out that the conditions at which the pyrolysis is performed can affect the heat demand (or supply) of the reaction, thus influencing the reaction mechanism. Apparently, the main reaction pathway in flash pyrolysis is endothermic, and a value of 129 kJ/kgvolatiles has been proposed for the heat of reaction,4 reflecting the latent heat of vaporization of the primary tars. Accordingly, a term accounting for this endothermicity must be included in the energy balance, as indicated in eq 4, resulting in the equation
cp
dTp 6 ) [h (T - Tp) + σ(Tem4 - Tp4)] dz dpFpvp gp g dLH2O dLDC evap ∆HH ∆HR (27) O 2 dz dz
which clearly shows that contributions from the latent heat of water evaporation and the heat of reaction are
Ind. Eng. Chem. Res., Vol. 41, No. 24, 2002 5997
Figure 9. Calculated and experimental weight loss at three temperature levels. Single first-order kinetics including the reaction heat (∆HR ) 129 kcal/kg), solid line. The prediction with ∆HR ) 0 of Figure 8 is also reported as a dashed line. The best-fit parameters are ln(A) ) 23.65 (A in s-1), E/R ) 16 128 K, L ∞dc ) 0.54.
regulated by the rate of water release and the rate of reaction, respectively. Equation 27 must, thus, also be solved together with eqs 18 and 23, which provide the rates just mentioned. The results are shown in Figure 9. To highlight the differences, Figure 9 reports also the predictions without the heat of reaction as a dashed line. Apparently, the predictions in this case are very close to those obtained neglecting the heat of reaction, which might explain why it is difficult to identify the magnitude and sign of the heat of the pyrolysis reaction using SFOR models. The SFOR model accounting for the heat of pyrolysis tends to improve the agreement with both short-residence-time data, which gives the lowest weight loss, and the results at higher yield, where the importance of the heat of reaction is more perceptible. Note that fitting the kinetics parameters partially compensates for the differences introduced by the heat of reaction. The most relevant conclusion, however, is that including the endothermicity of the decomposition reaction gives a best-fit value for the activation energy of 134 kJ/mol, which is very close to the literature value (140 kJ/mol) frequently found in previous studies at high heating rates that used SFOR models.4 We believe that this agreement demonstrates the consistency of our experimental investigation with previous laboratory studies, possibly involving different techniques, but we still note the limitations and evolution of the simple SFOR model. One of the more evident insufficiencies of the above prediction is the underestimation of the measured weight loss for short residence times (z ) 0.4 m), whereas in all other cases (z > 0.4), the model prediction fits in the experimental variance. Note that the heat of reaction taken from Milosavljevic and Suuberg4 is quite a significant figure, so that the particle temperature is markedly affected, remaining lower than that of the gas, as shown in Figure 10. Note, however, that the isothermal region, seen from the point of view of the particle transformation (i.e., conversion), is much wider than expected from the diagrams in terms of z. Perhaps the most significant limitations of SFOR models can be seen in W ∞DC (or L∞), which is determined as a constant parameter, even though it is known to depend on the experimental conditions. Consequently, the value identified here is of little use for scale-up of the process, even with the same cellulose. The same conclusion can be drawn for the kinetics parameters A
Figure 10. Temperature profiles of the gas (solid lines) and particles (symbols) with respect to conversion expressed as total weight loss, L. Single first-order kinetics including the reaction heat (∆HR ) 129 kcal/kg).
and E, though less intuitively, even if literature abounds with discussions of the different values identified by various researchers, all using the same SFOR model. Importance of the Activation Step The particular behavior exhibited by the experimental data that we obtained at the lower temperatures (400450 °C) can be well described with a two-step model for cellulose pyrolysis 1
DC 98 γV1 + (1 - γ)AC
(M1)
2
AC 98 βV2 + (1 - β)Ch The first activation step consists of a dramatic reduction of the degree of polymerization, yielding the socalled active cellulose. The value of γ is rather small, about 0.05, so the first step causes a limited weight loss, while the second step is the one that actually determines most of the weight loss. Note that we did not specify the nature of the volatiles produced, which can include light gases of condensable species (tars). According to the mechanism above, we can write the mass balances of the solid species as
dWDC ) -k1WDC dz
(28)
dWAC ) (1 - γ)k1WDC - k2WAC dz
(29)
dWCh ) (1 - β)k2WAC dz
(30)
vp vp
vp
where WAC and WCh are the masses of active cellulose and char, respectively. Defining the following yields, by analogy with eq 17
ξAC ≡
WAC W 0DC
ξCh ≡
WCh W 0DC
(31)
we can rewrite the material balances in eqs 28-30 in terms of the weight loss of dry cellulose and the yields of the two solid products
5998
Ind. Eng. Chem. Res., Vol. 41, No. 24, 2002
dLDC ) k1(1 - LDC) dz
(32)
dξAC ) (1 - γ)k1(1 - LDC) - k2ξAC dz
(33)
dξCh ) (1 - β)k2ξAC dz
(34)
vp vp
vp
includes an activation step, with limited weight loss, in addition to the loss of water initially contained in the cellulose. We next considered the existence of two parallel and independent reaction steps following the activation of cellulose according to the mechanism 1
DC 98 AC 2
The mass balances in eqs 32-34 are solved together with the energy balance of eq 27 and the rate of water loss expressed by eq 18. Once the functions LDC(z), LH2O(z), ξAC(z), and ξCh(z) are determined, the total weight loss L to be compared with the experimental data can be calculated as
L≡
AC 98 Ch 3
AC 98 V
In this case, only reaction 3 leads to significant weight loss. The material balances for all of the species involved are written as
W-W0 ) W0 0 (W H + W 0DC) - (WH2O + WDC + WAC + WCh) 2O
dWDC ) -k1WDC dz
(36)
dWAC ) 0.95k1WDC - (k2 + k3)WAC dz
(37)
vp vp
0 WH + W 0DC 2O
LDC + LH2O - ξAC - ξCh
(M2)
(35)
vp
dWCh ) k2WAC dz
(38)
The results of fitting all of the experimental data, from 400 to 600 °C, are shown in Figure 11. The twostep model contains a larger number of parameters because there are two Arrhenius constants, k1 and k2, but also because the stoichiometry (i.e., γ and β) was determined according to the data. The calculated value of γ is 4.6%, confirming that the first step causes a limited weight loss. The identification of a separate step is also crucial in explaining the experimental data at lower temperature (400 and 450 °C). The SFOR model cannot describe the data at 400 and 450 °C with the same set of parameters determined with the hightemperature data, which we saw above to be consistent with literature values. Moreover, SFOR predictions of our short-residence-time data at low temperature are rather poor, because of the steplike shape of the profiles that can be seen in Figure 11. Low-temperature measurements behave quite differently, as already observed,2 because the increase of conversion with residence time is very sharp for short residence times and then does not continue to increase significantly. Such behavior can be consistent with a mechanism that
vp
dWV ) k3WAC dz
(39)
)
1+
0 χH 2O
Two independent paths leading to volatiles and char clearly offer many more opportunities to explain a larger variety of data or, in other words, introduce more parameters that can be seen as additional degrees of freedom to fit the data. Calibration of the kinetic constants of the two reaction steps can be done precisely on the basis of experimental determinations of the amounts of char and volatiles produced. Unfortunately, in most experimental studies on cellulose pyrolysis only the loss of solid mass is measured. This is not a limitation in discriminating the roles of the two paths, because only the conversion of active cellulose to volatiles causes a loss of solid weight. This provides a constraint for the reactivity ratio of active cellulose to volatiles or char, i.e., k2/k3. We begin by simply observing that the ratio of the volatiles and char mass balances (eqs 38 and 39) yields
k3 WV ξV ) with ξV ≡ 0 ξCh k2 W
(40)
DC
because neither volatiles nor char is in the feed
W 0V ) 0
W 0Ch ) 0
Equation 40 constrains the kinetic parameters of the second and third steps, A2, A3, E2, and E3, to the volatiles and char yields as
ln
Figure 11. Calculated and experimental weight loss at all available temperatures. Two-step model with activation of cellulose. The best-fit parameters are ln(A1) ) 57.6 (A1 in s-1), E1/R ) 32 082 K, γ ) 0.046; ln(A2) ) 25.4 (A2 in s-1), E2/R ) 17 997 K, β ) 0.52.
k3 E3 - E2 1 ξg ) ln ) (ln A3 - ln A2) ξCh k2 R T
(41)
Neglecting the weight loss due to the activation step and the water loss, so that the function L can be determined by step 3 alone, and assuming that only char and volatiles leave the reactor, we obtain the following condition
Ind. Eng. Chem. Res., Vol. 41, No. 24, 2002 5999
ln
ξg ξg ) ln ) ξC 1 - ξg ln
(
)
W 0dc - WCh W 0DC L ) ln (42) 0 WCh 1-L W DC
Now, supposing that the entire process takes place in the isothermal part of the reactor, which is also the higher-temperature region in each experiment, from eqs 41 and 42, we would expect a plot of ln[L/(1 - L)] vs 1/T to be linear. The same analysis was performed on the data of Graham.5 The results are shown in Figure 12. Note that only the measurements obtained at z ) 1.2 m for the five temperatures tested were used because of the assumptions made, which require the maximum allowable length of the isothermal region. Excellent agreement between the two experimental investigations was found, closely matching the linearity expected. A single point of Graham,5 indicate by the × in Figure 12, behaves differently, but we cannot provide any explanation for this divergence from the information given by the authors.5 From the equation of the interpolating line in Figure 12, we can derive two conditions for the kinetic parameters of steps 2 and 3, namely
ln(A3) - ln(A2) ) 10.8
(43)
(E3 - E2)/R ) 9040 K
(44)
and
The rigorous solution of all of the energy and material balance equations, i.e., eqs 18, 27, and 36-39, and the fitting of the solution to the experimental data provides values for the kinetic parameters of each step, that is Aj and Ej, for j ) 1-3. The results are shown in Figure 13. Again, we observe that the activation step allows us to describe both the low- (400, 450 °C) and the high(500-600 °C) temperature data. The inclusion of the activation step and the particle temperature profile causes the parameters of the second and third steps to deviate from the constraints identified above, yielding
ln(A3) - ln(A2) ) 4.6
(45)
and
Figure 13. Calculated and experimental weight loss at all available temperatures. Two-independent-step model with activation of cellulose. The best-fit parameters are ln(A1) ) 58.9 (A1 in s-1), E1/R ) 32 862 K; ln(A2) ) 21.3 (A2 in s-1), E2/R ) 14 995 K; ln(A3) ) 25.9 (A3 in s-1), E3/R ) 18 861 K.
Figure 12 should use a convenient mean temperature. However, the influence of the heating period is the same for each residence time for a given T set, so that linearity is maintained. Interestingly, the extrapolation of the line in Figure 12 to low temperature (say 250 °C) predicts negligible weight losses. TGA tests can result in yields of volatiles higher than 50% at these low temperatures as well, but with residence times much larger than those obtained in the EFRs (approximately 0.5 s or less). There are at least two explanations for this discrepancy, both referring to more complex reaction mechanisms. One is that the short residence time allowed by the EFR precludes further reactions in series, e.g., interactions between volatiles and solid, that can increase the weight loss, as predicted by the multistep model of Broido et al.7 1
DC 98 AC 2
3
4
AC 98 D 98 E 98 char 5
AC 98 tar
(M3)
Note also that Broido et al.’s mechanism7 includes an activation step. A second possible explanation is a mechanism with one more parallel reactions following the activation steps. 1
(E3 - E2)/R ) 3866 K
DC 98 AC (46) 2
This difference arises mainly because the process taking place in our EFR is not exactly isothermal and
AC 98 char 3
AC 98 V1 4
AC 98 V2
Figure 12. Correlation between weight loss and temperature that provides additional constraints for the kinetic parameters. Experimental data labeled as this study refer to Visentin’s paper.2
(M4)
where the third and fourth steps can have very different time constants. The possibilities offered by a multistep mechanism are intriguing, but more detailed experimental data are required to validate them. We are presently working in this direction to obtain more details about the spectrum of products in our experimental setup. Note that more detailed kinetics mechanisms are expected to better approach the actual paths at the molecular level (intrinsic kinetics), thus enlarging the reliability of the predictions to a wider range of
6000
Ind. Eng. Chem. Res., Vol. 41, No. 24, 2002
Figure 14. Isothermal TGA (1 K/min, TMAX ) 250 °C) simulated through mechanism M1 and fitted with a single pseudo-first-order reaction (SFOR). Ea ) 234 kJ/mol.
Figure 15. Isothermal TGA (20 K/min, TMAX ) 350 °C) simulated through mechanism M1 and fitted with a single pseudo-first-order reaction (SFOR). Ea ) 115 kJ/mol.
conditions, effectively making the model a design tool for industrial-scale applications. Do Two Regimes Really Exist in Cellulose Pyrolysis? It has been observed4,28 that the many values of the activation energy in a pseudo-first-order model tend to cluster around two values, depending on the experimental conditions used. Experiments with high heating rates reaching a temperature above 600 K apparently result in low activation energies (∼140 kJ/mol), whereas experiments with slow heating (
