Estimation of gas phase mixing in packed beds

Estimation of gas phase mixing in packed beds

Combustion and Flame 153 (2008) 137–148 www.elsevier.com/locate/combustflame Estimation of gas phase mixing in packed beds S. Frigerio a , H. Thunman...

546KB Sizes 0 Downloads 34 Views

Combustion and Flame 153 (2008) 137–148 www.elsevier.com/locate/combustflame

Estimation of gas phase mixing in packed beds S. Frigerio a , H. Thunman b,∗ , B. Leckner b , S. Hermansson b a CMIC Dipartimento di Chimica, Materiali e Ingegneria Chimica, Politecnico di Milano, Piazza Leonardo da Vinci 32,

20133 Milan, Italy b Department of Energy Conversion, Chalmers University of Technology, SE-412 96 Göteborg, Sweden

Received 16 April 2007; received in revised form 11 May 2007; accepted 11 May 2007 Available online 12 July 2007

Abstract An improved model is presented for estimation of the mixing of gaseous species in a packed bed for fuel conversion. In particular, this work clarifies the main characteristics of mixing of volatiles and oxidizers in a burning bed of high-volatile solid fuel. Expressions are introduced to represent the active role of degradation of the solid particles in the mixing within the gas phase. During drying and devolatilization the solids modify the behavior of the gas flow: the volatiles released from the surface of the particles increase the turbulence in the system, and hence the rates of the homogeneous reactions under mixing-limited conditions. Numerical experiments are carried out to test the validity of this conclusion regarding mixing in different geometries. The flow of volatiles leaving the fuel particles is shown to contribute significantly to mixing, especially at low air flows through a bed. However, the fraction of the particle surface where volatiles are released and its orientation in the bed should be better determined in order to increase the accuracy of the estimates of turbulent mixing. © 2007 The Combustion Institute. Published by Elsevier Inc. All rights reserved. Keywords: Grate combustion; Fixed bed conversion; Packed bed conversion; Mixing; Fuel conversion

1. Introduction Packed-bed systems have been adopted for a wide range of applications. Recently, conversion of solid, high-volatile fuel, such as biomass, in fixed beds has become of interest, requiring a more detailed description of the process inside a bed of burning particles. Although differences obviously exist, some characteristics of such beds are similar to those of catalytic systems: in particular, both consist of a large number of randomly arranged immobile particles, passed * Corresponding author. Fax: +46 31 772 3592.

E-mail address: [email protected] (H. Thunman).

by a gas flow through the voids between the particles. Several models for catalytic fixed-bed systems [1–3] allow prediction of flow field, mass and heat transfer between gas and solid phase, and pressure drops across the bed. Even if the models of catalytic systems describe many phenomena in a burning bed, in noncatalytic systems additional expressions are required for drying, devolatilization, gasification, and combustion modeling. Such expressions take into account modifications of the solid (shrinkage, swelling, and fragmentation) and associated variations in the properties of the gas flow. The solids release volatile species (H2 O, H2 , CO, CO2 , light and heavy hydrocarbons) into the void space between the particles. The momentum of the volatile species pro-

0010-2180/$ – see front matter © 2007 The Combustion Institute. Published by Elsevier Inc. All rights reserved. doi:10.1016/j.combustflame.2007.05.006

138

S. Frigerio et al. / Combustion and Flame 153 (2008) 137–148

Nomenclature A a,b C Cmix D d Eact F f g h∗ k K0 L lc m n P Pe R Rdev Rkin Rmix Re S Sc t T uV U

Surface . . . . . . . . . . . . . . . . . . . . . . . . [m2 ] Pressure drops coefficients Concentration . . . . . . . . . . . . . [mol/m3 ] Mixing rate constant Gas dispersion coefficient . . . . . [m2 /s] Diameter . . . . . . . . . . . . . . . . . . . . . . . [m] Activation energy . . . . . . . . . . . . [J/mol] Flow rate. . . . . . . . . . . . . . . . . . [kg/m2 s] Mass flow . . . . . . . . . . . . . . . . . . . . [kg/s] Acceleration of gravity . . . . . . . . [m/s2 ] Height of active volume . . . . . . . . . . [m] Turbulent energy . . . . . . . . . . . . [m2 /s2 ] Pre-exponential factor . . . . . . . . . . [s−1 ] Pressure drop length . . . . . . . . . . . . . [m] Turbulent length scale . . . . . . . . . . . . [m] Mass . . . . . . . . . . . . . . . . . . . . . . . . . . . [kg] Number . . . . . . . . . . . . . . . . . . . . . . . . . [–] Gas pressure . . . . . . . . . . . . . . . . . . . . [Pa] Péclet number Pe = Uφ dP /D Gas constant . . . . . . . . . . . . . [J/(mol K)] Devolatilization rate . . . . . . . . . . . . [s−1 ] Arrhenius rate . . . . . . . . . . . . . . . . . [s−1 ] Mixing rate . . . . . . . . . . . . . . . . . . . [s−1 ] Reynolds number Re = ρU dP /μ Mass source term . . . . . [mol m−3 s−1 ] Schmidt number Time . . . . . . . . . . . . . . . . . . . . . . . . . . . . [s] Temperature. . . . . . . . . . . . . . . . . . . . . [K] Volatiles gas velocity at the particle surface . . . . . . . . . . . . . . . . . . . . . . . [m/s] Superficial gas velocity . . . . . . . . . [m/s]

vides energy, which is then dissipated within the void space, increasing the turbulence and thus the mixing rate. Consequently, the models for catalytic systems have to be extended with models for conversion of fuel in fixed beds, some of which deal with highvolatile fuels, e.g., [4–7]. A complete overview of the process includes both molecular and macroscale aspects; molecular diffusion, complemented by turbulence and mixing, should be taken into account when considering homogeneous reactions. The fuel (volatiles) and the oxidizer (air transported through the reactor) should be mixed in the reaction zone before reacting. The chemical kinetic rates are important for the ignition of the conversion process, but when the oxidative process has started, the exothermic reactions raise the temperature, and mixing can limit the conversion. The interaction between kinetics and mixing is often expressed by comparing the dominant mechanism among the homogeneous reactions and mixing. The rate is determined from the lowest

Uφ V V1,V2 x,y

Interstitial gas velocity . . . . . . . . . [m/s] Volume . . . . . . . . . . . . . . . . . . . . . . . . [m3 ] Volatiles inlets Coordinate directions

Greek letters β  ε ϑ λ μ ν ρ τ φ χ ϕ

Local pressure drop coefficient Difference Turbulent dissipation rate . . . . [m2 /s3 ] Fraction of the total particle surface located transverse to the fiber orientation Air-to-fuel ratio Gas viscosity . . . . . . . . . . . . . . . [kg/m s] Gas laminar viscosity . . . . . . . . . [m2 /s] Density . . . . . . . . . . . . . . . . . . . . [kg/m3 ] Tortuosity Void fraction of the bed Packing coefficient Particle shape factor

Subscripts F J L M m O P S T V

Fuel Junction Longitudinal Main flow Molecular Oxygen Particle Solid Turbulent conditions Volatiles flow

one of the conversion and mixing rates, S = min[Skin , Smix ],

(1)

where S is the source term in the mass transport equation. Skin is related to the kinetics and Smix = Rmix min{CO , CF } is the mixing rate Rmix times the minimum concentration of fuel or oxygen in the location of the reaction. In the present work the influence of the volatiles on the mixing is of interest. The significance of mixing in general was investigated by Magnussen and Hjertager [8]. They suggested that the mixing rate Rmix depends on the turbulence properties of the flow, both the turbulent kinetic energy k and the dissipation of turbulent energy ε, ε Rmix = Cmix , k

(2)

where Cmix is a mixing rate constant. This concept can be applied in a burning bed as well. A first attempt to estimate the turbulent properties in a (fluidized) bed

S. Frigerio et al. / Combustion and Flame 153 (2008) 137–148

for fuel conversion was made by Palchonok et al. [9] (following Deckwer [10], who worked with bubble columns). Palchonok et al. stated that “all the energy supplied into the furnace with the fluidization gas dissipates into heat due to turbulent motion of the phases.” Frictional energy loss is small and was neglected. They defined the dissipation rate related to the gas volume in the bed, ε ∝ ρP g(1 − φ)

P U U = , ρφ L ρφ

(3)

and the kinetic energy, k = (εlc )2/3 ,

(4)

where ρ is the density of particles (index p) and gas, φ the voidage of the bed, P the pressure drop across the bed height L, and U the superficial velocity. The length scale of turbulence is lc . The same concept was applied to a fixed bed by Yang et al. [11], using the Ergun equation [12] for the pressure drop P and relating the dissipation to a gas occupying the entire volume of the bed (or, alternatively, allowing the voidage to be included in the mixing rate coefficient). These authors treated mixing due to the main flow, but the additional impact of the volatiles is missing. This impact on the global mixing is certainly related to the stoichiometric condition of the process. In particular, during gas production (air-to-fuel ratio λ < 0.4), the effect of the volatiles’ momentum can be stronger, relatively seen, than in the case of combustion (λ > 1) because of the difference in total gas flow through the bed. Furthermore, the devolatilization rate depends on temperature and heating rate, affecting the instantaneous velocity of the volatiles. The objective of the present work is to derive a model to determine mixing, taking the effect of volatiles into account. The expression for mixing is validated by numerical simulation. The validity of the model depends on the velocity field: in a packed bed the deviation from strictly laminar flow becomes significant at a particle Reynolds number of Re > 10 [13], and therefore the velocity corresponding to Re = 10 is assumed to be the lower limit (Ulow ) of validity for the present model. To avoid fluidization, the upper limit is the minimum fluidization velocity Umf , and standard operating conditions of a bed for gasification or combustion of solid fuel are within the velocity range Ulow < U < Umf . This range is still valid when the Re number becomes smaller as a consequence of a rise in temperature.

2. Theory There are three aspects to be considered for modeling of the dispersion and the source terms in the

139

transport equation related to mixing and kinetics in a gas flow passing a burning bed of high-volatile fuels: • Molecular motion, depending on the species involved; • Turbulence in the main flow, the flow that transports the oxidizer from the inlet to the outlet of the bed, incorporating the volatile species released by the solids and later downstream also combustion intermediates and products; • Enhancement of turbulence due to the release of volatiles from the solids during drying (H2 O vapor) and devolatilization (light and heavy gaseous hydrocarbons, CO, H2 , CO2 , H2 O). In the following, the role of these items is analyzed. 2.1. Mixing due to the main flow The volatiles driven off from the hot particles should mix with the oxidant transported by the main flow before reacting. Many researchers have analyzed mixing effects on chemical reactions, starting with the work of Danckwerts [14]. The situation is analogous to the reaction between nonpremixed reactants, where the interaction between mixing and kinetics must be taken into account. In particular, a gas phase reaction is a process on a molecular scale and micromixing should be considered. On the other hand, under turbulent conditions the mixing on larger scales (mesomixing) can be controlling, and the micromixing can be neglected [15]. Mesomixing is a phenomenon that occurs on scales that are larger than those of molecular mixing but smaller than the reactor dimension, such as the size of the void regions between the particles [16]. The mixing is governed by the time scale for large eddies as described by the eddy-breakup model [17]: it depends on the ratio of the dissipation rate of turbulent energy ε and the turbulent kinetic energy k. The turbulence variables arise from the power transferred to the fluid, which is proportional to the ratio of the third power of the flow velocity, u, and a characteristic length, L, dependent on the geometry of the system. As the velocity depends on the mass flow through the bed, the volatiles released from the par ticles increase the mixing rate in two ways, as shown in Fig. 1: 1. The entering mass, transported by the main flow FM (the vertical arrow in Fig. 1), is increased by the additional mass flow of volatiles FV . The velocity of the main flow UM , and consequently the turbulence, rises through the bed. 2. FV brings energy depending on the third power of velocity, u3V , which becomes incorporated into

140

S. Frigerio et al. / Combustion and Flame 153 (2008) 137–148

Fig. 1. Schematic gas flow pattern (FM = mean flow rate, FV = flow rate of volatiles) inside a fixed bed of active and nonactive particles. The darker particles are not yet reacting.

the motion of the large eddies. uV is the velocity of the volatiles at the surface of the particle. Its value and direction depend on the arrangement of the particles and on the devolatilization rate. This additional effect from the volatiles on the generation of turbulence is usually neglected in the modeling of mixing. Spalding’s analysis [17] was used by Magnussen and Hjertager [8]. In their original work the constant Cmix = 4 was introduced into Eq. (2) to fit measurements. Later work [18], dealing with direct numerical simulation of homogeneous systems, indicates that a value of unity is more correct. For quantitative analysis this value is important, but for the purpose of the present work the distinction can be omitted, and a value of unity will be used in the following. The net rate of each homogeneous reaction is then determined as the lower of the kinetic and mixing rates (Eq. (1)). Assuming that Eq. (2) is valid for fixed-bed reactors, the source term is determined by Eq. (1). The application of the above relations requires the determination of the turbulence variables ε and k, using only such macroscopic properties of the bed as are available in standard models of fixed-bed reactors. The main procedure for solving this problem in connection with fixed-bed combustion [11] estimates the mixing in the main flow in terms of the macroscopic properties of the reactor. The pressure drop P in the vertical direction of a fixed bed is given by the sum of two terms: a viscous energy loss term, proportional to the fluid velocity, and an inertial loss (kinetic energy) term, proportional to the square of the velocity, 2. P /L = aUM + bUM

(5)

For porous beds, the coefficients a and b are commonly estimated by Ergun’s equation [12], P /L = 150

μ(1 − φ)2 ϕ 2 dP2 φ 3

UM + 1.75

ρ(1 − φ) 2 U , ϕdP φ 3 M (6)

Fig. 2. Simplified conduit geometry (right-hand side) representing the enlargements (with length L2 and width d2) and the restrictions (with length L1 and width d1) of the flow pathway of the gas flow among the particles in the fixed bed (left-hand side).

where ϕ is the shape factor of the particles and μ the gas viscosity. The diameter of a particle with an arbitrary shape is dP = 6VP /(AP ϕ). In the cases of spherical and square particles ϕ = 1. Based on Eqs. (2), (3), (4), and (6) the mixing rate becomes [11] Rmix,M = 150

(1 − φ)2/3 Dm lc φ lc

+ 1.75

(1 − φ)1/3 UM . lc φ

(7)

(The constant Cmix = 0.83 was introduced [11] to fit a set of simulations performed with a computational code.) The second term on the right-hand side (expressing the deviation from laminar flow) comes directly from Ergun’s equation, while the first one (the molecular term) is obtained in analogy with the struc = D /l ture of Eq. (6). A molecular velocity UM m c was introduced instead of the convective velocity UM , although it cannot be directly explained on the basis of the Ergun equation. Later, Thunman and Leckner [19] suggested that only the second term on the righthand side is related to mixing, since the first term in the Ergun equation is dominated by pressure drop due to friction, which is not transformed into turbulent movements, and a molecular contribution to mixing can be neglected in the turbulent regime. It is possible to demonstrate the validity of such an assertion by considering the work of Niven [13] on pressure drop inside a porous bed. In particular, Niven [13] proposed a simple model (Fig. 2) to describe the path inside the bed, showing that the gas mixing is due to flow separation behind each particle and to enlargements and restrictions, as well as to changes in the flow direction within the packed bed. The fixed bed is represented in Fig. 2 as a conduit of alternating restrictions and enlargements, whose dimensions are d1–L1 and d2–L2. The comparison between the simulated pressure drop of this config-

S. Frigerio et al. / Combustion and Flame 153 (2008) 137–148

uration and the one calculated from Ergun’s equation showed a satisfactory agreement. The pressure drops are determined by both viscous and form resistances (Eq. (5)), but the viscous one is converted into heat and therefore only form resistances affect the turbulence. From Niven’s analysis this contribution becomes ρ(1 − φ) 2 UM , P /L = χ dP φ 3

(8)

where χ depends on the tortuosity and packing arrangement, as well as on the shape of the particles. Accepting Eq. (6), which is valid for porous beds, the pressure drop can be written as P /L = 1.75

ρ(1 − φ) 2 UM dP φ 3

(9)

and the turbulent dissipation of energy becomes, according to Eq. (3), ε=

(1 − φ) 3 U P = 1.75 U . φρL dP φ 4 M

(10)

This relationship is independent of viscosity, since for macro structures viscous effects are negligible and the energetic quantities include only the inertial contribution. Using Eq. (4) for k and Eq. (10) for ε, the mixing rate of the main flow becomes Rmix,M =

ε1/3

ε = 2/3 (εlc )2/3 lc

(1 − φ)1/3 = 1.2 1/3 2/3 UM . dP lc φ 4/3

(11)

As a first assumption lc = dP . The new expression differs from the one proposed by Yang et al. [11], (Eq. (7)), because there is no viscous term, the constant is 1.2 instead of 1.75, and the term φ 4/3 appears instead of φ. 2.2. Mixing due to moisture and volatiles flow leaving particles: New expression Following the assumption that the power provided by the turbulent motion and then dissipated is the product of the pressure drop in the flow rate of the gas (Eq. (3)), the first step in taking care of the impact of volatiles consists in estimating the velocity of the volatiles leaving the particles’ surface. The devolatilization rate can be described locally inside a particle or for a small particle as a first-order reaction by an Arrhenius expression, dmV = −mV Rdev , (12) dt where Rdev = K0 exp(−Eact /RTS ) and mV is the remaining content of volatiles. In the case of a single

141

symmetric particle, the velocity depends on the particle’s surface AP and on the rate of devolatilization dmV /dt: uV =

dmV /dt . ρV AP

(13)

In a bed consisting predominantly of thermally large particles, drying and pyrolysis fronts move gradually into a particle, with temperature TS mostly controlled by the thermal conductivity of the particle. The process in the bed becomes smeared out: a conversion front with a certain thickness is observed in the bed. Because of the steep rate–temperature relationship (the exponent), drying and pyrolysis follow each other quite closely inside a particle [20]. Therefore, for the present purpose, the effective rates of drying and devolatilization can be deduced together, assuming a conversion front with thickness h∗ in the bed and using experimental data on conversion rate in the bed, FV , expressed as mass per square meter bed area Areact . The velocity of the volatiles averaged over the total surface of the devolatilizing particles in the bed, uV,bed , is then estimated as the ratio of the total mass of moisture and volatiles to the total surface in the active volume: uV,bed =

FV Areact . ρV AP,h∗

(14)

The total surface AP,h∗ of the devolatalizing particles is equal to that of a single particle AP times the number of particles n0P,h∗ in the volume Areact (1 − φ)h∗ ,  AP,h∗

= n0P,h∗ AP =

 Areact (1 − φ) ∗ h AP . VP

(15)

The flow field is influenced by the local velocity uV at the surface of a particle and not by the average velocity uV,bed . The local velocity depends on the fraction ϑ of the total surface of a fuel particle that is located transverse to its fibers, through which the gas is expected to be emitted, as illustrated in Fig. 3. When the fraction ϑ is smaller than unity, the velocity increases compared to the average velocity. This

Fig. 3. T-pipe schematization (right-hand side) of volatiles released horizontally from a square particle (left-hand side).

142

S. Frigerio et al. / Combustion and Flame 153 (2008) 137–148

In this case, the characteristic length, related to the pressure drop, is half of the total distance between the inlets of the volatiles L = x/2, assuming perfect symmetry between the two opposite inlets (see Fig. 3), while the characteristic length of turbulence is lc = x. The total pressure drop is twice that caused by each inlet (since the inlets provide the same energy). The associated mixing rate of the volatiles’ flow is obtained according to Eq. (2) with Cmix = 1,

(a)

Rmix,V = (b) Fig. 4. (a) Arrangement of two particles with adjacent surfaces releasing gas. Fiber orientation is (a) horizontal for both particles A and B; (b) horizontal for particle A and vertical for particle B.

implies that the volatiles are driven off from the particle mainly in the fiber direction and hardly across it. In case there is no predominant fiber direction, ϑ = 1. The release velocity of the volatiles into the bed becomes in general uV = uV,bed /ϑ.

(16)

Moreover, in a fuel bed the particles may be located in various ways and this irregular position of the particles also has to be taken into account. For instance, two-particle arrangements, either with adjacent releasing surfaces (Fig. 4a) or with one particle that releases gas against an inactive surface (Fig. 4b), can be strongly different with respect to total energy as well as its direction of release. These features are included in the factor ϑ , to be determined empirically, that expresses the amount of active surface area that is directed perpendicular to a gas channel passing the bed. When the volatile gases leave a particle they change abruptly their direction, causing a pressure drop that can be calculated like that of the flow in a channel at a T-pipe junction, as illustrated in Fig. 3, assuming a square particle and a horizontal direction of the flow of volatiles. The local pressure drop of a T-junction can be calculated by u2 PJ = β V ρV = u2V ρV , 2

PJ uV , L ρV

k = (εlc )2/3 .

(20)

2.2.1. Overlapping of contributions from the main flow and gas from the particles A semiquantitative description of the overlapping contributions to mixing from volatiles and main flow can be made, starting from the following assumptions: • The total mixing is the combination of the contributions from FM and FV . • The FM mixing is increased by the action of the FV . • The role of FV depends on the ratio of the momentum transported by the volatiles to that transported by the main flow. It is related to the respective, mass flows, applied to a specific area AF , calculated as fV ∝ ρV uV AF,V

and fM ∝ ρM UM AF,M . (21) AF depends on the particle arrangement. For the main flow (vertical direction) it corresponds to the empty area between the particles. The definition of the specific area for the volatiles is more critical and has to be determined empirically (being dependent on ϑ ). • The contribution of the volatiles α is then proportional to the ratio α=

fV . fM + fV

(22)

The resulting mixing rate obtains contributions from Eq. (11) (main flow, index M) and from Eq. (19) for the volatiles’ flow (index V): Rmix = Rmix,M + αRmix,V .

(23)

(17)

where the coefficient β is set equal to 2 [21]. From Eqs. (3) and (4) the turbulent variables ε and k are determined: ε=2

ε ε1/3 41/3 uV . = = 2/3 x (εlc )2/3 lc

(18) (19)

3. Evaluation 3.1. Strategy The resulting mixing rates can be evaluated simply by knowledge of the properties of the bed (x, dP , and φ) and with the respective velocities UM and uV :

S. Frigerio et al. / Combustion and Flame 153 (2008) 137–148

143

Fig. 5. Two-dimensional (A, B, C, D, E) and three-dimensional (F) bed arrangements of uniform spheres (A, B, C, F) and squares (D, E).

by Eq. (11) for the main flow through the bed and by Eq. (20) for the mixing rate produced by the volatiles emitted in the reaction zone. These expressions need validation. In the absence of direct measurements, the validity of the expressions can be assessed by comparing with Eq. (2), where the quantities k and ε are obtained by a numerical model of the corresponding gas flow through a packed bed. Finally, the significance of the contribution of the volatiles will be estimated using the situation in an actual biomass fuel bed undergoing conversion, ranging from devolatilization to combustion. However, first, the consistency of the numerical model will be verified by comparing modeling results from various configurations of two- or three-dimensional beds (2D, 3D) with available measurements of dispersion of gas in packed beds. 3.2. The numerical model Fluent 6.1.2 is employed for the numerical solution of the steady state continuity, momentum, and species balance equations, including the standard k–ε model equations. The default parameters were used as defined in the Fluent manual [22]. Only the Schmidt number (default value 0.7) is modified (see below). An upwind first-order scheme was used in all cases. The pressure correction scheme followed the SIMPLE algorithm. The most important assumptions are as follows: • 2D and 3D geometries with either spherical or square particles. • The flow is considered stationary, even if transient phenomena can affect the process. • For the present purpose it is not necessary to treat combustion directly. Instead combustion is assumed to produce a certain temperature. Hence, isothermal conditions are assumed in the computational domain. • The fluid is incompressible. • A near-wall model describes the flow near the particles’ surface; the near-wall region, affected by viscosity, is resolved all the way down to the viscous sublayer.

• In all simulations Re > 100. The definition of mixing requires that both k and ε be predicted under well-developed turbulent conditions. • The standard k–ε and RSM turbulence models were compared, but the differences were small, and only results from the former are reported here. 3.3. The configurations simulated The numerical code was applied on several packed bed configurations, illustrated in Fig. 5. Configurations A, B, and C are 2D beds of circular particles with various packings, whereas D and E are similar cases with square particles. The case F is a 3D bed of spherical particles. An analysis, analogous to those made for the beds, is performed to estimate turbulence directly related to the flow of volatiles from two particles. The model, shown in Fig. 6, is similar to that used by Niven [13] (Fig. 2). The release of volatiles is simulated by two couples of mass-flow inlets (V1 and V2) that differ in the distances (x1 = 2x2 ). Only one couple of inlets (two particles) is used at a time. This configuration allows straightforward and independent estimation of the most important features of the mixing: • mixing due to the main flow (excluding the volatiles’ outflow);

Fig. 6. Simplified channel geometry, representing the gas flow between square particles in the fixed bed.

144

S. Frigerio et al. / Combustion and Flame 153 (2008) 137–148

(a)

(b)

Fig. 7. (a) Rate of drying and devolatilization (FV ) related to dry mass of the fuel versus air flow (Fair ) [23], expressed as mass rate per square meter bed area. Line: stoichiometric conditions. (b) Maximum temperature in the combustion front [23]. In both part open symbols are cylindrical particles of length 34 mm and diameter 8 mm [!], 12 mm [e]; filled symbols [", a] (A–F) are measurements used in the present work as test cases.

(a)

(b) Fig. 8. Packed bed. Ratio of the derived (Eq. (11)) to the simulated (Eq. (2)) expression for ε/k (a) for spherical particles with geometries A, B, C, F (Fig. 5) and (b) for square particles ε/k with geometries D, E (Fig. 5). Main flow.

• mixing due to the flow of volatiles (excluding the main flow as well as the volatiles entering V2 when considering V1, and analogously for V1 when considering V2). Finally, a practical case provides data for evaluation of the significance of release of gas from the particles in applications: a fixed bed with 0.5 × 0.5 m cross section, operated in counter flow (the combustion front moves against the air flow), and fired with wood of various sizes (8 and 12 mm) [23]. The measured weight of this reactor gives the total release of gases from the particles. The flow of volatiles leaving the fuel particles can be calculated according to Eqs. (14)–(16), and relevant temperatures for evaluation of physical properties inside a reaction front are

Fig. 9. Single conduit. Ratio of the derived (Eq. (11)) to the simulated expression (Eq. (2)) for ε/k in main flow in zones V1 and V2 (Fig. 6).

shown in Fig. 7: the rate of drying and devolatilization is determined from the mass loss of fuel (Fig. 7a) and the temperature (Fig. 7b) versus air flow (expressed per square meter bed area) for substoichiometric, stoichiometric, and superstoichiometric conditions. The stoichiometric line in Fig. 7a represents the amount of fuel converted if all the oxygen in the air is consumed. The diagram shows that in the substoichiometric range (left-hand side), more fuel is converted than that corresponding to the oxygen supplied: moisture and volatiles from the solid fuel leave the bed together with the combustion products. Owing to substoichiometry, these volatiles are only partly burned in the mixing region. Only a small part of the char is burned or gasified, while the major part accumulates above the reaction front that proceeds downward in this batch-fired combustor. In the combustion zone (right-hand side), all fuel is burned, and the excess oxygen leaves with the combustion gases. The computational grids for the various configurations mentioned were chosen to accurately represent the wall regions in all configurations and then expanded in the channel spaces in such a way that a reduction in grid size would not further influence the results. This means that the number of computational cells varied from case to case from 8000 to 30,000 for the cases in Fig. 5 and from 6000 (main flow) to 200,000 (volatiles flow) in Figs. 6, 8, and 9.

S. Frigerio et al. / Combustion and Flame 153 (2008) 137–148

145

The input flow was chosen to be even over the cross section of the bed. The boundaries of particles as well as bounding walls were treated as smooth walls with a no-slip condition. The particular flow patterns at the side walls, such as in the configurations of Fig. 5, were not considered specifically, as the turbulence data of interest were chosen in a representative channel inside the beds far from the side walls.

4. Validation of the dispersion model Dispersion of gas in a packed bed is a classical topic of investigation and a great deal of measurement data are available. This information is useful in the present context to verify the ability of the numerical model to reproduce the flow through the packed bed. It is particularly important to compare the dispersion coefficients resulting from the simulation with those from the classical measurements of the longitudinal Péclet number [24,25] represented by 1 1 1 1 = + , PeL τ Pem 2

(24)

with a tortuosity τ = 1/φ 0.5 related to bed voidage φ. The value of 2 in the right-side term is the asymptotic PeL reached for Re > 5 or Pem > 10 and hence is representative in the present context. Cases A to F were used in the comparison. The numerical code considers dispersion in the form of molecular and turbulent terms [22]. The latter term is determined as DT =

μT . ρScT

(25)

Turbulent viscosity is the independent variable modeled by the code instead of turbulent diffusivity, and the turbulent Schmidt number (ScT = νT /DL ) has to be known. Using the default setting ScT = 0.7 leads to overestimation of the asymptote 2 for the longitudinal Péclet number in all cases tested. This fact could imply that the simulated dispersion is lower than what would be achieved from measurement (PeL = 2). However, the adequacy of the adopted default value ScT = 0.7 to predict turbulence in a flow between particles has not been experimentally assessed. In fact, ScT has been determined in several investigations and is found to vary in the range 0.2–1, depending on the geometry and on the distance from the wall (Hinze [26] obtained 0.625 for the core region of a turbulent pipe-flow; He et al. [27] suggested adopting 0.2 for a jet in crossflow; Koeltzsch [28] reported a polynomial expression to fit experimental data within the boundary layer and 0.4 for the outer region). In the present work, the adoption of an ScT of about 0.4

Fig. 10. Comparison between simulated data (symbols) and expression (Eq. (24), curve) for the longitudinal Péclet number, for geometries A–F, Fig. 5. Also, the overestimated range of value obtained using the default value ScT = 0.7 is reported (dotted region).

gives a satisfactory agreement with the experimental data, as seen in Fig. 10. This value will be used in the following evaluation. In fact, this is a calibration of the dispersion model and implies that the turbulent viscosity in Sc has been adjusted to fit the average relation (Eq. (24)) of a great quantity of independent measurements of dispersion in packed beds.

5. Results and discussion 5.1. Packed bed: Mixing due to the main flow The analysis of the main flow in packed beds, represented by the configurations in Fig. 5, consists of a comparison between the mixing calculated by the numerical simulation and that of the mixing model. The simulated values are mass-weighted averages of the ratio ε/k according to Eq. (2), while the calculated value is given by Eq. (11). The results are shown in Fig. 8 as a function of Reynolds number for the flow through the bed. First, Fig. 8 shows that in all cases the calculated and simulated mixing rates are of the same order of magnitude (as reported by Yang et al. [11]). Second, the agreement is more evident in the case of strong turbulence, that is, either for Re numbers well above 100 or using square particles. This can be related to the ability of the code to reproduce highly turbulent flows. There is no obvious difference between the 2D and the 3D cases, but the 3D case happens to agree very well with the calculated value. The conclusion is that 2D simulations give acceptable results, and they are used in all simulations except Case F in Fig. 8a. 5.2. Single conduit 5.2.1. Mixing due to main flow The second simulation investigates the mixing in the simplified configuration of Fig. 6. First, only mix-

146

S. Frigerio et al. / Combustion and Flame 153 (2008) 137–148

Fig. 11. Single conduit. Ratio of the derived (Eq. (20)) to the simulated expression (Eq. (2)) for ε/k in zones V1 and V2 due to volatiles only (Fig. 6).

ing due to the main flow is investigated. Both inlets for volatiles in Fig. 6 are simulated as walls, and the main flow is characterized by the vertical component (y-direction). The mixing rate is then evaluated with the same analysis as above for the porous bed. The results from the numerical code and Eq. (2) are compared with the ones from Eq. (11) in Fig. 9, considering the constriction between each couple of inlets as the area of interest. The comparison confirms the good agreement shown in the previous analysis for both V1 and V2 (in particular V2, characterized by a stronger turbulence). 5.2.2. Mixing due to volatiles only The validity of Eq. (20) for the estimation of mixing due to volatiles is tested with the geometry of Fig. 7. An estimate of the velocity of volatiles is made using Eq. (12) with the kinetic values reported in [20]. These values are appropriate to preserve the range of operability of the turbulence models of the numerical calculations. The comparison between calculated and simulated mixing is quite satisfactory, as shown in Fig. 11. Since the volatiles flow out from the surfaces, there is no direct relation between the Reynolds number of the vertical flow in a channel through the bed (as typically defined) and that of the flow of volatiles. Turbulence intensity (fluctuating velocity component compared to average velocity of the gas) is therefore reported instead of Re on the horizontal axis. This comparison shows that, despite using a simple model (T-pipe analogy), the mixing rate, expressed by the ratio ε/k, is well reproduced. Adoption of a characteristic length scale is quite optional, even if its value most likely is close to the one chosen (the largest dimension of the void space). In this work lc is always the largest dimension between the particles. 5.3. Qualitative mixing estimates Now, as the ability to estimate the mixing due to main and volatile’ flow is demonstrated, a qualitative comparison between the two contributions can be made for the situation in real fuel-bed converters

(ranging from gasifiers to combustors) with representative data taken from the example in Fig. 7. The reaction front can be assumed to be on the order of a few particle diameters. Three diameters is chosen as a representative number, as shown by measurements carried out simultaneously with those presented in Fig. 7 [23], and in agreement with calculated estimates [19]. This gives the height (h∗ = 3dP ) of the volume of interest (to be used in Eq. (15)). The measured air flows are contained between Ulow and Umf . In the most critical case (smallest particle size), dP = 0.008 m, Ulow ∼ = 0.05, and Umf = 0.7 m/s. It is now possible to demonstrate the magnitudes of mixing due to volatiles and to main flow. Three conditions of flow are considered for both dp = 0.008 m (points A, B, and C) and 0.012 m (points D, E, F) in Fig. 7). 5.3.1. Volatiles mixing The mixing due to volatiles only is evaluated by Eq. (20), adopting ϑ = 0.5 when calculating uV by Eq. (16). Density is calculated assuming a molar mass of 0.028 kg/mol (corresponding to CO) and an average temperature between the maximum temperature reported in [23] (Fig. 7b) and 600 K, adopted as the lower temperature required for devolatilization [29]. It is evident from Fig. 7b that temperatures are higher (1500 K) in the case of superstoichiometric (combustion) conditions than during gasification (1200–1300 K). 5.3.2. Main flow mixing The main flow velocity is simply calculated from the air flow entering the reactor at the same temperature as the volatiles UM =

Fair . ρair φ

(26)

The mixing rate is then evaluated according to Eq. (11). The results are reported in Fig. 12 for the points A to F indicated in Fig. 7. The global mixing is higher for smaller particles, according to Eqs. (11) and (20). Moreover, the contribution of the volatiles to mixing is of the same order of magnitude as that of the main flow and should be taken into account in combustion calculations, at least in the low and intermediate velocity ranges. In particular, in the case of the lowest main flow rate (gasification conditions, points A and D), mixing due to volatiles is a significant part of the total mixing. When the main flow rate is higher, the two effects are still comparable, although under combustion conditions, the contribution of the main flow is clearly dominant. The final goal of the analysis is a quantitative estimation of global mixing. At the moment, the lack of information on the empirical parameter (ϑ ) prevents

S. Frigerio et al. / Combustion and Flame 153 (2008) 137–148

147

Fig. 12. Comparison between mixing due to main ( ) and volatile’ flow (2). Cases A, B, C and D, E, F from Fig. 7 refer to diameters of 8 mm and 12 mm, respectively. Data are taken from Fig. 7.

such an analysis, but the effect of the volatiles has been clearly demonstrated.

6. Conclusion A formulation that includes the effect of volatiles is proposed for mixing in the gas phase in packed fuel beds. Furthermore, an improved model for the estimation of the reaction in the gas phase in a fuel bed is presented: the reactions cannot be estimated by kinetic analysis only (as for well-mixed gases), but the possible contribution of mixing has to be considered. Moreover, the model recognizes that the mixing is caused by energy provided both by the main flow and by the volatile species released from the particles during the solid degradation (drying, devolatilization, gasification, and combustion processes). In particular, expressions are provided to estimate the effect of both the main flow and the volatiles’ outflow on the mixing rate. An acceptable agreement was found, as verified by comparing the estimates from the expressions suggested with results determined by simulation of the flow with a commercial numerical code for various geometries. The comparison with and without volatiles shows the significance of taking the outflow of volatiles into account in the determination of mixing of gas species during analysis of conversion of solid fuel in a fixed bed, especially under substoichiometric conditions where the main gas flow is small. Although the results are ready for application in packed bed analyses, the direction of the release of volatiles is not yet well described, owing to the anisotropy of the biomass material and the disordered orientation of the particles in a bed. Therefore, a better estimate of these features of bed conversion is still needed to improve the accuracy of the model.

References [1] G.F. Froment, K.B. Bischoff, Chemical Reactor Analysis and Design, second ed., Wiley, New York, 1990. [2] V. Balakotaiah, S.M.S. Dommeti, Chem. Eng. Sci. 54 (1999) 1621–1638.

[3] H.P.A. Calis, J. Nijenhuis, B.C. Paikert, F.M. Dautzenberg, C.M. van den Bleek, Chem. Eng. Sci. 56 (2001) 1713–1720. [4] C. Di Blasi, Chem. Eng. Sci. 55 (2000) 2931–2944. [5] D. Shin, S. Choi, Combust. Flame 121 (2001) 167–180. [6] Y.B. Yang, H. Yamauchi, V. Nasserzadeh, J. Swithenbank, Fuel 82 (2003) 2205–2221. [7] H. Thunman, B. Leckner, Fuel 82 (2003) 275–283. [8] B.F. Magnussen, B.H. Hjertager, Proc. Combust. Inst. 16 (1976) 719–729. [9] G.I. Palchonok, F. Johnsson, B. Leckner, in: M. Kwauk, J. Li (Eds.), Circulating fluidized bed technology, in: Proc. of the 5th Int. Conf. on Circulating Fluidized Beds, Science Press, Beijing, 1996, pp. 440–445. [10] W.-D. Deckwer, Chem. Eng. Sci. 35 (1980) 1341– 1346. [11] Y.B. Yang, Y.R. Goh, V. Nasserzadeh, J. Swithenbank, Mathematical Modelling of solid combustion on a travelling bed and the effect of channelling, presented of the Third Int. Symp. on incineration and flue gas treatment technologies, Brussels, 2001. [12] S. Ergun, Chem. Eng. Prog. 48 (1952) 9–94. [13] R.K. Niven, Chem. Eng. Sci. 57 (2002) 527–534. [14] P.V. Danckwerts, Chem. Eng. Sci. 8 (1958) 93. [15] S. Valerio, M. Vanni, A.A. Barresi, Chem. Eng. Sci. 49 (1994) 5159–5174. [16] J. Baldyga, J.R. Bourne, S.J. Hearn, Chem. Eng. Sci. 52 (1997) 457–466. [17] D.B.M. Spalding, Proc. Combust. Inst. 13 (1970) 649– 657. [18] V. Eswaran, S.B. Pope, Phys. Fluids 31 (1988) 506– 520. [19] H. Thunman, B. Leckner, Proc. Combust. Inst. 30 (2005) 2939–2946. [20] H. Thunman, K. Davidsson, B. Leckner, Combust. Flame 137 (2004) 242–250. [21] A. Cocchi, Elementi di termofisica generale ed applicata [Elements of general and applied thermophysics], Soc. Ed. Esculapio, Bologna, 1990. [22] Fluent 6.1 User’s Guide. [23] M. Axell, Förbränningsförlopp i en bädd av biobränsle (Processes of combustion in a bed of biofuels), Thesis for Licentiate of Engineering, Department of Energy Conversion, Chalmers University of Technology, Göteborg, Sweden, 2000. [24] D.J. Gunn, Trans. Amer. Inst. Chem. Eng. 47 (1969) 351–359. [25] R.H. Wilhelm, Pure Appl. Chem. 5 (1962) 403–421.

148

S. Frigerio et al. / Combustion and Flame 153 (2008) 137–148

[26] J.O. Hinze, Turbulence, McGraw–Hill, New York, 1975. [27] G. He, Y. Guo, T. Hsu, The effect of Schmidt number on turbulent scalar mixing in a jet-in-crossflow, Int. J. Heat Mass Transfer 42 (1999) 3727–3738.

[28] K. Koeltzsch, Atmos. Environ. 34 (2000) 1147– 1151. [29] E. Biagini, L. Tognotti, Proc. of Joint Meeting of Scandinavian-Nordic and Italian Sections of the Combustion Institute, Ischia, 2003.