3852
Ind. Eng. Chem. Res. 2005, 44, 3852-3861
Steam-Gasification of Coal in a Fluidized-Bed/Packed-Bed Reactor Exposed to Concentrated Thermal RadiationsModeling and Experimental Validation Peter von Zedtwitz† and Aldo Steinfeld*,†,‡ Department of Mechanical and Process Engineering, ETHsSwiss Federal Institute of Technology Zurich, ETHsZentrum, CH-8092 Zurich, Switzerland, and Solar Technology Laboratory, Paul Scherrer Institute, CH-5232 Villigen, Switzerland
The steam-gasification of coal in a fluidized-bed or a packed-bed chemical reactor is considered using an external source of concentrated thermal radiation for high-temperature process heat. The energy equation that couples heat transfer with the chemical kinetics is solved by means of a numerical model that incorporates the Monte Carlo ray-tracing technique for nonisothermal, nongray, absorbing, emitting, and scattering media. The reaction kinetics are described by Langmuir-Hinshelwood type rate laws. Validation is accomplished by comparing the numerically computed temperature profiles, product gas composition, and reaction extent with the experimentally measured values using a tubular quartz reactor directly exposed to high-flux irradiation. For the packed bed, the temperature increases monotonically because the internal radiative exchange approaches a conduction-like heat transfer within the bed. For the fluidized bed, the temperature increases rapidly in the first one-quarter of the bed and then reaches a constant value because of the strong fluidization in the upper bed region derived from the 5-fold volumetric growth due to gas formation and thermal expansion. Above 1450 K, the product composition consisted mainly of an equimolar mixture of H2 and CO, a syngas quality that is notably superior than that typically obtained in autothermal gasification reactors (with internal combustion of coal for process heat), besides the additional benefit of the upgraded calorific value. 1. Introduction Hybrid solar/fossil endothermic processes, in which fossil fuels are used exclusively as the chemical source for H2 production and concentrated solar power is used exclusively as the energy source of process heat, offer a viable route for fossil fuel decarbonization.1,2 An important example of such hybridization is the endothermic steam-gasification of coal to synthesis gas (syngas). The syngas product, besides being a high-quality fluid fuel for efficient combined-cycles and fuel cells, is cleaner than its solid feedstock because its energy content has been upgraded by the solar input in an amount equal to the enthalpy change of the reaction. Hybrid solar/ coal processes, such as the one presented in this paper, create a transition path toward solar hydrogen.2 Several chemical aspects of the gasification of solid carbonaceous materials are summarized in the literature.3-7 In conventional gasifiers, the energy required for heating the reactants and for the heat of the reaction is supplied by burning a significant portion of the feedstock, either directly by internal combustion or indirectly by external combustion.7 Internal combustion, as applied in autothermal reactors, results in the contamination of the gaseous products, while external combustion, as applied in allothermal reactors, results in lower thermal efficiency because of the irreversibilities associated with indirect heat transfer. Alternatively, the advantages of supplying solar energy for process heat are three-fold: (1) The calorific value of the feedstock * Corresponding author. E-mail:
[email protected]. † ETHsSwiss Federal Institute of Technology Zurich. ‡ Paul Scherrer Institute.
is upgraded. (2) The gaseous products are not contaminated by the byproducts of combustion. (3) The discharge of pollutants to the environment is avoided. The gasification of carbonaceous materials and related reactions have been performed using concentrated solar energy in exploratory, early studies with coal, oil shales, biomass, and other carbon-containing feedstocks.8-15 Several solar reactor concepts have been proposed and tested with small-scale prototypes.16-18 More recently, the reaction kinetics of the steamgasification of coal were investigated for a fluidized-bed reactor directly exposed to an external source of concentrated thermal radiation.19 The direct irradiation of particle suspensions provided an efficient means of heat transfer directly to the reaction site, as corroborated as well in previous experimental studies.20-23 Numerous studies have been performed in the area of numerical modeling of fluidized-bed systems. The modeled thermochemical processes included gasification, pyrolysis, and combustion of carbonaceous materials (coal, biomass, emulsions, etc.).24-27 In particular, the net-flow modeling from the emulsion to the gas phase for bubbling and pressurized fluidized-bed gasifiers allowed for better prediction of experimental results.28,29 For processes carried out at >1300 K, radiative transfer within the fluidized particulate media and among the fluidized bed and reactor walls becomes a predominant mode of heat transfer.30-34 This paper presents the modeling of a directly irradiated reactor, its experimental validation, and its use to predict the reactor’s performance under various operational conditions. The modeled chemical reactor consists of a quartz tube containing a fluidized bed or packed bed of coal particles undergoing steam gasification. An
10.1021/ie050138w CCC: $30.25 © 2005 American Chemical Society Published on Web 04/28/2005
Ind. Eng. Chem. Res., Vol. 44, No. 11, 2005 3853
Figure 1. Geometrical configuration: (a) ETH’s high-flux solar simulator with tubular quartz reactor; (b) elements of tubular quartz reactor and fluidized bed.
external source of concentrated thermal radiation delivers the required high-temperature process heat. The energy conservation equation is formulated for the gasparticle system that couples radiation, conduction, and convection heat transfer to the chemical reaction. The Monte Carlo (MC) ray-tracing method is applied to solve the 3D radiative exchange within the particulate bed and within the reactor walls. This powerful numerical technique enables the treatment of problems beyond the reach of analytical techniques and is flexible enough to remove diffuse-gray-isothermal-isotropic simplifying assumptions for treating nondiffuse, nongray, nonisothermal, and anisotropic problems involving complex 3D geometries.35 Examples of recent applications include the modeling of pulverized coal fired furnaces,36 real gas mixtures,37 particle clouds having spectral and directional dependent properties,38-40 and directly irradiated solid-gas systems undergoing thermochemical transformations.41,42 Previous modeling studies for fluidized coal gasifiers include the net-flow concept.
Figure 2. Spectral intensity of emitted radiation from the Ar arc.
origin at the center of the trough elliptical mirror, as shown in Figure 1. 3. MC Analysis
2. Geometrical Configuration The geometrical configuration is shown schematically in Figure 1. The radiation source is the Swiss Federal Institute of Technology’s (ETH’s) high-flux solar simulator.43 This research facility provides a rapid external source of intense thermal radiation in the visible and IR spectrum that approaches the heat transfer characteristics of highly concentrating solar systems. It consists of a high-pressure Ar arc enclosed in a 27-mm diameter, 200-mm length water-cooled quartz envelope. The arc is close-coupled to a precision elliptical trough reflector, with one of the linear foci coinciding with the arc. The focal plane of the solar simulator is, thus, defined as the horizontal plane containing the second linear focus. The elliptical mirror is truncated 9.2 cm above the focal plane to permit external access. With this arrangement, radiative power fluxes exceeding 4500 kW/m2 are attained at the focal plane and confined within a 45° rim angle. Power, power fluxes, and temperatures can be adjusted to meet the specific requirements by simply varying the electrical input power to the arc electrodes. A tubular quartz reactor containing a packed bed or fluidized bed of coal particles is positioned vertically, with half its height at the focal plane. The system of coordinates is selected with the
Radiation emitted by the Ar arc is assumed to be isotropic. Thus, the direction of a generic ray is given by the unit vector uˆ ) [sin θ cos φ, sin θ sin φ cos θ]T, randomly determined by the cone angle and circumferential angle (θ and φ, respectively)
φ ) 2πRφ
(1)
-1
θ ) cos (1 - 2Rθ) where R is a random number chosen from a uniform set [0, 1]. Two spectral distributions of the emitted radiation are employed: Planck’s distribution for a 5780 K blackbody, to simulate the solar spectrum, and the spectral distribution for the Ar arc, shown in Figure 2.43 For Planck’s distribution, the associated wavelength is found from
λ ) F0-λT-1(Rλ)/T
(2)
for T ) 5780 K. For the Ar arc, the following equation is used
Rλ )
∫0λiλdλ/∫0∞iλdλ
(3)
3854
Ind. Eng. Chem. Res., Vol. 44, No. 11, 2005
absorbed, a count is recorded at the location of absorption, and its history is terminated. To ensure no energy accumulation in radiative equilibrium, a new ray is isotropically emitted from the same location at a wavelength found by solving the implicit equation
∫0λaλeλb(λ,T)dλ Rλ ) ∞ ∫0 aλeλb(λ,T)dλ Figure 3. Spectral absorption coefficient of a quartz layer.
with iλ given by Figure 2. Rays undergo none, single, or multiple reflections at the mirror elements. The direction of the reflected ray rˆ is calculated using a combined model for diffuse/specular reflection that incorporates an angular error for mirror imperfections.44 For diffuse reflection, the direction is chosen randomly from a set that is weighted according to Lambert’s cosine law. For specular reflection, the direction is calculated by
rˆ ) uˆ - 2(nˆ uˆ )nˆ
(4)
where nˆ is the normal unit vector of the surface at the point of reflection, nˆ ) ∇S/|∇S|. S(x,y,z) ) 0 is the equation of the reflecting surface composed of the arc, elliptical, lateral flat, and angular flat mirrors (see Figure 1). When a ray exiting the reflector hits the quartz walls of the tubular reactor, it is either reflected or refracted. Total reflection occurs when the angle between the incoming ray and the surface normal, θ1, is greater than θ1,max ) sin-1(n2/n1), where n1 and n2 are the refractive indices of the medium it exits and enters, respectively. Otherwise, if θ1 < θ1,max, the probability of reflection is given by Fresnel’s equation for the directional-hemispherical specular reflectivity
[
]
2 cos2(θ1 + θ2) 1 sin (θ1 - θ2) 1+ F′λ(λ,θ1) ) 2 sin2(θ + θ ) cos2(θ - θ ) 1
2
1
2
(5)
where θ2 is the angle between the refracted ray and the normal to the surface, calculated by Snell’s law:
θ2 ) sin-1
(
)
n1 sin θ1 n2
(6)
If F′λ > RF, the ray is reflected. In both cases of reflection, the direction of the reflected ray is determined using eq 4. If the ray is not reflected, it will be refracted in the direction
(
rˆ ) cos θ1
)
n1 n1 - cos θ2 nˆ + uˆ n2 n2
(7)
The quartz tube is modeled as a 3D nonisothermal, nongray, absorbing, emitting, nonscattering participating medium. The path length lq to the probable point of attenuation inside the quartz is determined by
lq ) -
lnRl aλ,q
(8)
where aλ,q is the quartz absorption coefficient, bandapproximated using measured data45 as shown in Figure 3. If lq is smaller than the distance the ray needs to travel until reaching the boundary, the ray is
(9)
where eλb(λ,T) is evaluated at the temperature of the location of emission. Equation 9 is solved by applying the band approximation for the three spectral intervals of Figure 3: [0, λ1], [λ1, λ2], and [λ2, ∞] with λ1 ) 0.17 µm and λ2 ) 3.0 µm, yielding
{
λ)
( )
RλaP a1 a1F0-λ1T for Rλ e aP RλaP - (a1 - a2)F0-λ1T F0-λT-1 a2 a1F0-λ1T (a1 - a2)F0-λ1T + a2F0-λ2T for < Rλ e (10) aP aP RλaP - (a1 - a2)F0-λ1T - (a2 - a3)F0-λ2T F0-λT-1 a3 (a1 - a2)F0-λ1T + a2F0-λ2T for < Rλ aP
F0-λT-1
(
)
(
)
where aP is the Planck mean absorption coefficient
aP(T) ) a1F0-λ1T + a2(Fλ1T-λ2T) + a3(1 - F0-λ2T)
(11)
and a1, a2, and a3 are the values of the absorption coefficients in the three spectral intervals, respectively. Once the direction and wavelength of the new emitted ray are known, the ray is tracked to determine the new location of attenuation. Note that if lq is greater than the distance to the boundary, the ray is either reflected or refracted at the quartz walls according to the procedure described earlier. A ray leaving the quartz layer to the outside is assumed lost. A ray leaving the quartz layer to the inside enters the fluidized/packedbed region [xbb, xbt] or the lower/upper void regions. The fluidized/packed bed is modeled as a 3D nonisothermal, nongray, absorbing, emitting, and anisotropically scattering participating medium. Each particle is assumed to be opaque, isothermal (a good approximation for Bi , 1), and spherical (a generally good assumption for most irregularly shaped, randomly oriented particles) and to have independent scattering (as justified by referring to the independent/dependent scattering regime map for the range of particle volume fractions and size parameters used in this study46). Radiative properties are computed using geometrical optics for large particles, since for an average particle diameter of 1 mm and for the spectral interval of interest, the size parameter ξ ranges between 300 and 3000. For simplification, the gas phase is taken to be a nonparticipating medium; its contribution to the radiative transfer has been estimated to be 5), diffusely reflecting spheres are35
φr ) 2πRφr
σλ ) Fλrp2πNp aλ ) (1 - Fλ)rp2πNp Φ(θ) )
(12)
8 (sin θ - θ cos θ) 3π
where Np is the particle number density and θ is the angle between the directions of the incoming and scattered rays. Fλ is the spectral hemispherical reflectivity of coal, obtained from measurements47 and shown in Figure 4. Note that the extinction coefficient in the bed region, κb ) σλ + aλ ) rp2πNp, is independent of the wavelength because, for large particles, it is linked to the free path length that the ray travels in the bed without interaction. Thus, the path length lb to the probable point of attenuation inside the bed is calculated from
lb ) -ln(Rl)/κb
(13)
At this location, another random choice, depending on the albedo, determines whether the ray is either absorbed or scattered. The probability of absorption is 1 - Ωλ, evaluated at the ray’s wavelength. If the ray is scattered, its direction forms an angle θ with the incident direction, obtained by numerically solving
Rθ )
If the ray is absorbed, a count is recorded at the location of absorption and its history is terminated. The absorbed energy is either transferred by convection to the gas stream or consumed by the endothermic gasification reaction. A portion of the absorbed energy is reradiated. To ensure there is no energy accumulation in radiative equilibrium, a new ray is emitted from the same location. The direction of isotropic emission is randomly determined by the polar and circumferential angles (θ and φ), given in eq 1. The new ray is emitted by the nongray medium at a wavelength found by solving eq 9, using the absorption coefficient from eq 12 and evaluated at the Tp of the location of emission. Once the direction and wavelength of the new emitted ray or of the old scattered one are known, the ray is tracked to determine the new location of attenuation. In approach 2, the fluidized bed is modeled as a set of randomly positioned spherical particles with a Gaussian particle size distribution. The MC technique follows a generic ray that interacts with the ensemble of particles. When a hit occurs, a random choice determines whether the ray is either absorbed or reflected at the particle’s surface, depending on the surface spectral absorptivity Rλ ) (1 - Fλ) evaluated at the ray’s wavelength. If the ray is reflected, the direction of diffuse reflection is chosen randomly from a set that is weighted according to Lambert’s cosine law,
2 ∫0θΦ(θ,ψ) sin θ dθ ) 3π [θ - 43 sin(2θ) +
1 2
1 θ cos(2θ) (14) 2
]
using the Newton algorithm. The scattering phase function Φ gives, in statistical terms, the probability of scattering in a certain direction. The azimuthal angle ψ of the scattered ray is randomly chosen between 0 and 2π, because Φ is independent of ψ for spherical particles. Since the scattering is assumed elastic, the wavelength of the scattered ray remains unchanged.
(15)
θr ) sin-1xRθr
where φr and θr are the polar and azimuth angles relative to the normal to the surface at the point of incidence. Otherwise, the procedure is equivalent to the one applied for the continuous media model. Eventually, each ray originated in the radiation source (Ar arc) that is incident in the reactor, or each ray emitted by the quartz medium or by the fluidizedbed medium, would undergo scattering/absorption by the participating media and reflection/refraction at the quartz boundaries. The history of a generic ray is then a complete random sequence, which terminates when it is absorbed by either the fluidized bed or the quartz walls or lost to the surroundings. This procedure is repeated for a large number of rays to obtain statistically meaningful results. 4. Chemical Kinetics The kinetic modeling is based on the analysis of elementary reaction mechanisms describing reversible adsorption/desorption processes and irreversible surface chemistry.19 Assuming plug flow conditions and pseudofirst-order reactions, a set of simplified rate laws of the Langmuir-Hinshelwood type with neglected inhibition terms are formulated to describe the formation and consumption of each gas species:
rH2O ) - K1pH2O rH2 ) K1pH2O rCO ) K1pH2O + 2K2pCO2 - 2K1K3pH2OpCO rCO2 ) - K2pCO2 + K1K3pH2OpCO
}
(16)
The reaction rate rj of gas species j (j ) H2, H2O, CO,
3856
Ind. Eng. Chem. Res., Vol. 44, No. 11, 2005
Table 1. Activation Energy and Frequency Factor for Arrhenius Relation (Eq 18) K1 K2 K3
EA
K0
205.3 kJ‚mol-1 9.01 kJ‚mol-1 -363.4 kJ‚mol-1
1220 mol‚kg-1‚s-1‚Pa-1 5.33‚10-7 mol‚kg-1‚s-1‚Pa-1 1.34‚10-19 Pa-1
Parc/Nrays. For the pseudocontinuous MC model (approach 1),
Q˙ iem,p ) 4aiP,pσTip4 dVp
(25)
For the particle-discrete MC model (approach 2),
and CO2) for heterogeneous surface reactions is defined as
Q˙ iem,p ) AipgRiP,pσTip4
1 dnj 1 dnj rj ) ) aS mC dt AS dt
where Aipg is the total particle surface in element i and aP and RP are the Planck mean absorption coefficient and the Planck mean absorptivity, respectively:
(17)
where mC is the mass of coal, AS is the area of the active surface, and aS is the specific active surface area. Note that the units of the rate constants are: for K1, mol‚s-1‚g-1‚Pa-1; for K2, mol‚s-1‚g-1‚Pa-1; and for K3, Pa-1. The system of four coupled differential equations (eq 16) is solved numerically by Runge-Kutta and by iterating on the values of K1, K2, and K3 to minimize the difference between the theoretically calculated and experimentally measured molar flow rates of the products. The temperature dependence of each Ki is determined by the Arrhenius relation:
( )
Ki ) K0,i exp -
EA,i RTp
(18)
The apparent activation energy and frequency factor obtained by linear regression are shown in Table 1. Since the Ki represent complex reaction mechanisms rather than elementary steps, a negative value for the apparent activation energy, as in the case of K3, is possible. 5. Energy Conservation For the purpose of discretization, the quartz tube is divided into a large number of rings, each assumed isothermal, and the particle suspension in the bed is divided into a large number of disks, each assumed isothermal. These subdivisions are illustrated in Figure 1 and are denoted with the superscript i (i ) 1...n). The boundary conditions are the temperature and molar 0 , n˘ 0Ar, n˘ 0CO composition of the incoming gas flow: T0g, n˘ H 2O 0 0 ) n˘ CO2 ) n˘ H2 ) 0 and total pressure ptot ) ¢. Steadystate energy conservation applied to each subsystem yields:
+
Q˙ icd,p
+
Q˙ imix
(19)
for gas: Q˙ ichm + Q˙ ipg + Q˙ iqg ) H˙ iout - H˙ iin
(20)
i ) Q˙ iem,q + Q˙ iqg + Q˙ iqa + Q˙ icd,q for quartz: Q˙ abs,q
(21)
The radiation terms are found by MC: i Q˙ abs,p ) wNip
(22)
i ) wNiq Q˙ abs,q
(23)
∫0∞ aλeλb(λ,T) dλ
1 σT4
∫0
1 RP ≡ σT4
∞
(27) Rλeλb(λ,T) dλ
The coefficient aP is computed using eq 11 for the quartz and using eq 9 for the bed. The conduction terms are:
for the quartz rings: Q˙ icd,q ) -Aqkq(Ti+1 - 2Tiq + Ti-1 q q )/dx (28) - 2Tip + Ti-1 for the bed: Q˙ icd,p ) -Abskp(Ti+1 p p )/ dx (for i ) 2...n - 1) Q˙ icd,p ) -Abskp(T2p - T1p)/dx (for i ) 1) Q˙ icd,p
)
Abskp(Tnp
+
Tn-1 p )/dx
(29)
(for i ) n)
with Ab and Aq being the bed’s and the quartz’s crosssectional areas. The particle mixing term, only considered for the fluidized-bed reactor, is given by i-1 i Q˙ imix ) n˘ C,mix (hC,Tpi - hC,Tpi-1) + n˘ C,mix (hC,Tpi - hC,Tpi+1) (30) i is the molar exchange of C between where n˘ C,mix elements i and i + 1 in both directions and hC,T is the specific enthalpy of C at temperature T. The convection terms are:
Q˙ ipg ) Aipg hipg(Tip - (Tig + Ti-1 g )/2)
(31)
Q˙ iqg
(32)
)
Aqghiqg(Tiq
Q˙ iqa
i ) Q˙ iem,p + Q˙ ipg + Q˙ ichm + for particles: Q˙ abs,p
Q˙ iC
aP ≡
(26)
)
-
(Tig
Aqahiqa(Tiq
+
Ti-1 g )/2)
- Ta)
(33)
Aqg and Aqa are the areas of the quartz-gas and quartzair contact surfaces, respectively. The heat transfer coefficients are calculated using the following correlations:
between particles and gas in a fluidized bed:48 Nupg ) 0.03Rep1.3 (34) between particles and gas in a packed bed:49 Nupg ) 2 + 1.1Prg1/3Rep0.6 (35)
(24)
between quartz and gas/particles in a fluidized bed:49 Nuqg ) 0.525Rep0.75 (36)
The variable w is the power carried by a single ray in an MC run with a sample of Nrays rays and w )
between quartz and gas/particles in a packed bed:49 Nuqg ) 2.26Rep0.8Prg0.33 e-6dp/dbdp/db (37)
Q˙ iem,q
)
i 4aP,q σT iq4
dVq
Ind. Eng. Chem. Res., Vol. 44, No. 11, 2005 3857
between quartz and surrounding air,50 evaluated at the film temperature: 2 0.387Raqa1/6 Nuqa ) 0.825 + (38) (1 + (0.492/Pra)9/16)8/27
[
]
Finally, the chemistry terms are
Q˙ ichm )
∑j n˘ r,ji hj,T
(39)
i p
i )(hC,T ip - hC,Ta) Q˙ iC ) -(n˘ ir,CO + n˘ r,CO 2
(40)
where hj,T is the specific enthalpy of species j at i is temperature T (j ) H2, H2O, CO, and CO2) and n˘ r,j i the mole number increase due to the reaction, n˘ r,j ) mi rij, with mi being the coal mass in element i. The kinetic rates are given by eq 16. The terms for the enthalpy flow in and out of the control volume are
∑j n˘ ij hj,T
i
g
-
∑j n˘ i-1 j hj,T
i-1 g
(41)
+ with the mole number of the outgoing flow n˘ ij ) n˘ i-1 j i . Note that the mole flow n˘ 0j entering the disk i ) 1 n˘ r,j is given by the boundary conditions. The system of nonlinear eqs 19-21 for one element was solved by the Newton algorithm, using ∆T ) 1 K and f ) 10-5 W. Five iterations were needed to converge with a difference of
