Sorption-Enhanced Steam

Article pubs.acs.org/IECR

Simulations of Steam Methane Reforming/Sorption-Enhanced Steam Methane Reforming Bubbling Fluidized Bed Reactors by a Dynamic One-Dimensional Two-Fluid Model: Implementation Issues and Model Validation Jannike Solsvik,* Rafael A. Sánchez, Zhongxi Chao, and Hugo A. Jakobsen* Department of Chemical Engineering, Norwegian University of Science and Technology (NTNU), N-7491 Trondheim, Norway ABSTRACT: Multiphase flows and in particular dispersed flows are frequently encountered in a variety of important process facilities in the chemical industry; herein, the fluidized bed reactors. Whereas the internal and external solid fluxes in the fluidized bed reactors are generally determined by empirical correlations or prescribed in the Kunii−Levenspiel type of models, the twofluid models serve as highly relevant in the progress of commercial reactors for processes such as the novel sorption-enhanced steam methane reforming (SE-SMR) technology because the solid flux incorporated in the two-fluid model allows for dynamic modeling of interconnected fluidized bed reactors. Nevertheless, the two- and three-dimensional two-fluid models are yet too computationally demanding for chemical process evaluation studies due to the complexity of the gas−solid flow in the bed. Hence, these models are not efficient for studies of the processes such as the SE-SMR technology where the solid particles are transported between reactor units for utilization and recovery of the characteristic solid property, i.e. CO2-capture and CO2release. Thus, the aim of the present study is to derive a model that allows for a more complex description of the fluidized bed reactors relative to the frequently used Kunii−Levenspiel type of models. On the other hand, the model should not predict details in the flow, as the two- and three-dimensional two-fluid models, in order to ensure reasonable simulation costs. Hence, in this study, the classical SIMPLE algorithm is extended to compressible two-phase reactive flows. The governing equations are cross-sectional averaged to smooth out the details in the flow. The suggested dynamic one-dimensional Eulerian−Eulerian twofluid model is applied to investigate the reactive gas−solid flows of the SMR and SE-SMR processes; and moreover, the simulation results are compared with the results of the more complex two-dimensional model with a solid stress closure based on kinetic theory of granular flow. The one-dimensional model predictions of the chemical process performance are in good agreement with the corresponding profiles predicted with the two-dimensional model. The deviations are larger comparing the internal flow details but these do not owe significant impact on the chemical process which to a large extent is determined by the imposed temperature in the reactor.

■

INTRODUCTION Steam methane reforming (SMR) is currently the predominant industrial route for hydrogen production. The development of alternative concepts for production of hydrogen via SMR has attracted considerable attention. A novel concept is the sorptionenhanced steam methane reforming (SE-SMR) process, which involves the addition of a solid sorbent into the SMR reaction system for the selective removal of CO2; and thereby, shifting the equilibrium toward increased hydrogen production. A regeneration step where the CO2 is released from the sorbent, and thus the capture activity is reintroduced, is necessary to make the SESMR process economically viable. Hence, the characteristic sorbent reactant utilized in the SE-SMR processes must be exposed to different reactor operation conditions in a cyclic manner. The SE-SMR process thus consists of two main steps: (i) reforming and (ii) regeneration. Because of the heterogeneous nature of the SE-SMR reactions and the relatively low chemical activity lifetime of the sorbent, the circulating fluidized bed is a suitable reactor concept that provides possible continuous operation of the process where the solid flux circulating between the reformer and regenerator reactor units aims to successively recover and utilize the chemical activity of the sorbent. © 2013 American Chemical Society

Although the SE-SMR technology is not commercialized, the concept of combining reformer reactions and separation in hydrogen production is not new.1 Rostrup-Nielsen2 reports that the first description of the conversion of hydrocarbons in the presence of steam and a CO2 acceptor was published in 1868. In 1933, Williams3 was issued a patent for a process in which steam and methane react in the presence of a mixture of lime and catalyst to produce hydrogen. However, there are many challenges associated to the commercialization of the SE-SMR technology. A sorbent regeneration step is necessary to make the SE-SMR process economically viable. The sorbent reversibility, i.e. decline of sorbent capacity over multiple carbonation/ calcination cycles, is a key factor affecting the process economy.4 Moreover, the appropriate sorbent should also fulfill a number of other properties like high sorption capacity, fast sorption/ desorption kinetic rates, and sufficient mechanical strength. Thus, efforts have been placed on investigation and development of, e.g., sorbent material,5−8 kinetic rate models for sorption,5,9 Received: Revised: Accepted: Published: 4202

December 4, 2012 January 27, 2013 February 19, 2013 February 19, 2013 dx.doi.org/10.1021/ie303348r | Ind. Eng. Chem. Res. 2013, 52, 4202−4220

Industrial & Engineering Chemistry Research

Article

product layer,10−13 and multiple carbonation/calcination cycles.5,13−16 On the other hand, relative to the conventional SMR process, the advantages of using a sorbent in the reforming process are that higher methane and CO conversions can be attained at lower temperatures and a product that gives higher H2 yield is possible. Additional benefits derivable from the SE-SMR process are the minimization of the coking potential, the elimination of the downstream H2 purification steps, the diminution of excess steam in the reforming operation, and the reduction of CO in the gas.17 The SE-SMR process operated in fluidized bed reactors has been described with different levels of complexity ranging from the simpler Kunii−Levenspiel type of models18,19 to rigorous multifluid models based on kinetic theory of granular flows (KTGF).20−23 The Kunii−Levenspiel type of models24−27 are frequently used in simulations of fluidized beds. This relatively simple model provides a complete formulation for solid distribution and mass balance. Adopting a Kunii−Levenspiel model for bubbling fluidized beds, the solids are commonly assumed to be stagnant in the reactor.24,28,29 Applying the Kunii−Levenspiel model for riser simulations, the solid distribution is prescribed in terms of empirical correlations instead of the solid phase continuity and momentum equations. Considering the circulating fluidized bed design, an intertransfered solid flux between the reactor units is required. Employing the Kunii−Levenspiel model for circulating fluidized beds, the solid flux is described or determined with addition correlations.30 The SE-SMR process is dynamic in nature as the solid density and thus the flow behavior like the solid flux changes with time and reaction performance. Thus, a hypothesis is that the Kunii− Levenspiel models are not appropriate to describe dynamic processes such as the SE-SMR technology. On the other hand, the two-fluid Eulerian−Eulerian models31 treat the gas and solid particles as interpenetrating continuous fluids. Hence, a solid flux is incorporated in the two-fluid model which allows for dynamic modeling of interconnected fluidized bed reactors. Lun et al.32 and Gidaspow33 have derived KTGF models which have been widely adopted for modeling and simulation of fluidized bed reactors, e.g. the cold-flow studies of Lindborg et al.34 and Wang et al.35 Two- and three-dimensional two-fluid models are yet too computationally demanding due to the complexity of the gas− solid flow in the fluidized bed reactors. Moreover, chemical reactions and the handling of large geometries of commercial reactors challenge the presently available computing capacity. In particular, the computations becomes excessively time-consuming for the additional consideration of interconnected reactors in a circulating fluidized bed design because the dynamic solid fluxes exchanging between the fluidized bed reactors need to be incorporated in the numerical solution algorithm. Such simulations are highly relevant and important in the progress of commercial fluidized bed reactors for processes such as the novel SE-SMR technology. A one-dimensional two-fluid model has the advantages of considerable lower computational costs relative to the two- and three-dimensional models. Whereas details in the complex gas−solid flow within a fluidized bed reactor; as obtained with the two- and three-dimensional models, is lost with a one-dimensional two-fluid model, the crosssectional averaged two-fluid model presents an improvement over the Kunii−Levenspiel type of models regarding the important solid fluxes transferred between the reactor units in a interconnected fluidized bed design.

The present study of the one-dimensional two-fluid model is performed to elucidate whether the model can be reasonably adopted for further simulations of interconnected fluidized bed reactors with a dynamic solid flux transferred between the reactor units. Dynamic solid circulation between units that operate at different conditions (e.g., temperature and feed composition) is an inherent requirement for the novel SE-SMR technology operated in fluidized bed reactors. A less computationally demanding one-dimensional model to study the performance of interconnected reactor units will be an important contribution to the progress of the commercialization of circulating fluidized bed reactors intended for the SE-SMR technology. The SE-SMR process can be operated using either a combined catalyst/sorbent pellet design unifying the catalytic and capture properties in a single pellet36−38 or by separating the catalytic and capture properties into two different pellet types. For simulation studies of the SE-SMR process using the two-pellet design, a multifluid Eulerian model may be derived giving the possibility for a catalyst phase, sorbent phase, and gas phase; e.g., the cold flow model of Chao et al.39 and the SE-SMR study by Carlo et al.40 and Chao et al.41 Recently, Solsvik and Jakobsen13 and Rout et al.42 performed a numerical study of a pellet holding both catalytic and capture properties with the objective to investigate the SE-SMR process. With this pellet design, the Eulerian model holds a single solid phase that interpenetrates with the gas phase. In the present study, the one-pellet design is adopted.

■

CHEMICAL REACTIONS Steam reforming of natural gas is the predominant production route to hydrogen for large-scale industrial applications. The kinetic reaction rates suggested by Xu and Froment43 are based on the following reactions: CH4(g) + H 2O(g) = CO(g) + 3H 2(g)

(1)

CO(g) + H 2O(g) = CO2 (g) + H 2(g)

(2)

CH4(g) + 2H 2O(g) = CO2 (g) + 4H 2(g)

(3)

In the SE-SMR process, the CO2 capture reaction by a CaObased sorbent is presented as CaO(s) + CO2 (g) = CaCO3(s)

(4)

The rate equation adopted in the simulations for the CO2 adsorption is taken from the work of Sun et al.9

■

ONE-DIMENSIONAL TWO-FLUID MODEL In the present study, a one-dimensional two-fluid model describing gas−solid flows with chemical reactions in fluidized Table 1. Constitutive Equations for Stresses Viscous stress ∂v τk = εgμg k (103) ∂z Wall friction stress

τs,wall =

fkf ρk vkvk 2

(104)

bed reactors is derived. In the continuum approach, the problem is formulated in terms of mass, species, heat, and momentum balances for each of the phases in an Eulerian reference frame. 4203

dx.doi.org/10.1021/ie303348r | Ind. Eng. Chem. Res. 2013, 52, 4202−4220


Article

Table 2. Constitutive Equations for Internal Heat Transfer Effective conductivity ⎡ μkt ⎤ ⎥Ck kkeff = ⎢kkm + ⎢⎣ ρk Pr t ⎥⎦

(105)

Molecular conductivity

kgstat

kgm =

εg kgstat

ksm =

εs

(1 −

εs )

(106)

(ϕA + (1 − ϕ Λ))

(107)

⎤ ⎡ A−1 B 2 A B−1 1 Λ= − (B + 1)⎥ ln − ⎢ 1 − B /A ⎣ (1 − B /A)2 A B 1 − B /A 2 ⎦

A=

ksstat kgstat

(108)

(109)

n

kgstat =

∑ yj ((j + )jT + *jT 2 + +jT 3) j

⎛ε ⎞ B = 1.25⎜⎜ s ⎟⎟ ⎝ εg ⎠

(110)

10/9

(111) −3

ϕ = 7.26 × 10

(112)

Continuity equation for phase k (= g,s)

Table 3. Constitutive Equations for Internal Mass Transfer Effective gas phase diffusivity ⎡ μgt ⎤ ⎥C Dg,effi = ⎢Dg,mi + Diff ⎢⎣ ρg Sc t ⎥⎦ (113)

∂ ∂ (εkρk ) + (εkρk vk) = Γk ∂t ∂z

56

Molecular diffusion coefficient (Wilke ) 1 − ωg, i 1 Dg,mi = ω Mg ∑ r = 1 n g, r MD r≠i

r ir

(114)

Momentum equation for gas phase

Binary diffusion coefficient57

Dir =

(5)

0.00266 × 101325Tg 3/2 0.012 pΩir Σir 2Mir1/2

(115)

∂ ∂ (εgρg vg) + (εgρg vgvg) ∂t ∂z ∂pg 4εgτg,wall ∂τg + εgρg gz + 4 g = −εg − + dt ∂z ∂z

The governing equations describing the reactive flow are presented in the sequent.

(6)

Table 4. Constitutive Equations for Interfacial Heat Transfer Interfacial heat transfer 6ε Q gsi = s hgs(Ts − Tg) dp

(116)

Interfacial heat transfer coefficient (spherical particles) kg hgs = [(7 − 10εg + 5εg 2)(1 + 0.7Rep0.2Pr1/3) + (1.33 − 2.4εg + 1.2εg 2) dp

Rep0.7Pr1/3]

(117)

Particle Reynolds number εgρg |vs − vg|d p Rep = μg (118) Prandtl number μg Cpg Pr = kg (119) 4204



Article

Table 5. Constitutive Equations for Interfacial Force Transfer Interfacial force

4 g = − 4 s = β(vs − vg)

(120)

The following constitutive equation for the interphase drag coefficient, β, is suggested by Gidaspow.33

Rep =

w=

(1 − εs)ρg |ug − us|d p μgm

(121)

a tan(150(εg − 1)) π

+ 0.5

(122)

⎧ 24(1 + 0.15Re 0.687) p ⎪ Rep < 1000 ⎪ Rep CD = ⎨ ⎪ ⎪ 0.44 Rep ≥ 1000 ⎩

(123)

⎛ C |U − U |ρ ε ε −2.65 ⎞ ⎛ ε 2μ m ρg |Ug − Us|εs ⎞ s g s g s g ⎟ ⎟ + w⎜ 3 D g β = (1 − w)⎜⎜a +b 2 ⎟ ⎜4 ⎟ d d d ε p p ⎠ ⎝ g p ⎝ ⎠

(124)

with a = 150 and b = 1.75. The weight w defined by eq 122 denotes an interpolation between the Ergun equation and the Wen58 drag closure.

In order to close the two-fluid model, constitutive equations are required for: (i) stresses (Table 1), (ii) internal heat transfer (Table 2), (iii) internal mass transfer (Table 3), (iv) interfacial heat transfer (Table 4), (v) interfacial momentum transfer (Table 5), and (vi) solid phase collision pressure (Table 6). To solve the mathematical model, the finite volume discretization technique is chosen.

Table 6. Constitutive Equations for Solid Phase Collision Pressure Solid phase collision pressure in terms of the modulus of elasticity: ∂ε ∂ (εsp ) = G(εg) s (125) ∂z s,coll ∂z with G(εg) given as

G(εg) = 10−10.46εg + 6.557

■

(126)

FINITE VOLUME METHOD The two-fluid model equations comprise a set of partial differential equations which require special techniques for

Momentum equation for solid phase ∂ ∂ (εsρs vs) + (εsρs vsvs) ∂t ∂z ∂pg 4εsτs,wall ∂τ ∂ε = −εs − + s + εsρs gz − G(εg ) s ∂z dt ∂z ∂z + 4s

(7)

Species mass balances for phase k (= g,s) ∂ ∂ (εkρk ωk , i) + (εkρk vkωk , i) ∂t ∂z ∂ωk , i ⎞ ∂⎛ eff = ⎜εkρk Dk , i ⎟ + Mi R k , i ∂z ⎝ ∂z ⎠

Figure 1. Grid for a one-dimensional problem. (8)

finding the solution. The idea is to convert the continuous equations into their discrete counterparts and end up with a system of algebraic equations that can be solved using standard equation solvers. In the finite volume method the partial differential equations are integrated over a domain subdivided into a number of cell volumes.31,44−46 Volume integrals of convective and diffusive terms in the partial differential equations are converted into surface integrals and evaluated as fluxes at the surface of each cell. The finite volume method is conservative by construction. Another property of the method is that numerous schemes and procedures can be design in order to solve the two-fluid model equations. On the other hand, even though the finite volume has become a very popular method deriving discretizations of partial differential equations, the method suffers from low accuracy and low convergence rates. By introducing a generalized variable φ that may represent, e.g., velocity, temperature, or species fraction, the conservative

Temperature equation for gas phase εgρg C p

∂Tg

+ εgρg C p vg

∂Tg

∂z ∂t ⎛ ⎞ T ∂ ∂ g )R rxSMR = + Q gsi ⎜εgkgeff ⎟ + ∑ ( −ΔHrxSMR i i ∂z ⎝ ∂z ⎠ rxi g

g

(9)

Temperature equation for solid phase ∂Ts ∂T + εsρs C pvs s s ∂t ∂z ∂ ⎛ eff ∂Ts ⎞ ⎜εsks ⎟ + ( −ΔH Cap)RCap − Q gsi = ∂z ⎝ ∂z ⎠

εsρs C p

s

(10)

It is noticed that the reformer reactions are allocated to the gas phase whereas the CO2-capture reaction take place in the solid phase. 4205



Article

components. Moreover, the momentum equation contains a contribution from the pressure, which has no analog in the generic equation. This term is treated as a surface force and thus approximated using the same quadrature rule as employed for the surface integrals in the generic transport equation. However, due to the close connection of the pressure and the continuity equation, the treatment of this term and the arrangement of variables on the grid play an important role in order to construct stable and accurate solution methods for the pressure−velocity coupling. Two-Fluid Model Equations. The basic discretization of the two-fluid model equations is similar to the approximations of the corresponding transport equations for single phase flow. A minor difference is that the two-fluid model equations contain the phase fraction variables that have to be approximated in an appropriate manner. More important, to design robust, stable, and accurate solution procedures with appropriate convergence properties for the two-fluid model equations, emphasis must be placed on the treatment of the interface transfer terms in the phase momentum, heat, and mass transport equations. Because of these extra terms, the coupling between the different equations is even more severe for multiphase flows than for single phase flows.31 Spalding47−51 suggested to calculate the phase fractions from the respective phase continuity equations in an approximate implicit volume fraction−velocity−pressure correction procedure. However, experience did show that it was difficult to conserve mass simultaneously for both phases with this algorithm. For this reason, Spalding suggested that the volume fraction of the dispersed phase may rather be calculated from a discrete equation that is derived from a combination of the two continuity equations. An alternative form of the latter volume fraction equation, particularly designed for fluids with large density differences, was later proposed by Carver.52 In this method, the continuity equations for each phase where normalized by a reference mass density to balance the weight of the error for each phase.31 Two techniques to handle the interphase coupling terms are to solve the momentum and heat balance equations of both phases simultaneously with a coupled solver, or, alternatively, the partial elimination algorithm (PEA) proposed by Spalding47,48 can be adopted; either applied with a segregated solver, or as in this study, in combination with a coupled solver.53 The working principle of the PEA is to weaken the strong coupling between the phases through a partial elimination of the variables in the interphase coupling term. The mathematical impact of the partial elimination algorithm on the solver for the resulting system of algebraic equations is to reduce the effect of the strong coupling between the phase momentum equations on the algebraic equations matrix condition number and thereby improve the convergence properties of the problem.31 The interfacial transfer fluxes are generally modeled as the product of the difference in the phase values of the primitive variables (driving force) multiplied by the transfer (proportionality) coefficients. A generic flux Ik can be expressed on the form:

differential equation used in the one-dimensional two-fluid formulation can be expressed by the following generic equation: ∂φk ⎞ ⎫ ∂ ∂ ∂⎛ (εkρk φk ) + (εkρk vkφk ) = ⎜εk Λk ⎟⎪ ∂t ∂z ∂z ⎝ ∂z ⎠ ⎬ ⎪ + S(φk ) + K (φl − φk ) ⎭ for l ≠ k

(11)

where constant cross-sectional area is assumed. Moreover, Λk is a generalized diffusion coefficient, Sk is a generalized source, and K is a generalized proportionality coefficient. Equation 11 is integrated over time and axial dimension:

∫Δt ∂∂t (∫Δz (εkρk φk) dz) dt + ∫ ∫ (εkρk vkφk ) dz dt Δt z ⎛

∂φ ⎞

∫Δt ∫z ⎜⎝εk Λk ∂zk ⎟⎠ dz dt + ∫ ∫ S(φk ) dz dt + ∫ ∫ Δt z Δt z =

dz dt

⎫ ⎪ ⎪ ⎪ ⎪ ⎬ ⎪ ⎪ K (φl − φk )⎪ ⎪ ⎭

for l ≠ k

(12)

Further, adopting the conventional notation; as illustrated in the one-dimensional staggered grid in Figure 1, the latter equation takes the discrete form: ⎫ ⎪ ⎪ ⎪ [εk n + 1ρk n + 1vk n + 1φk n + 1]e ⎪ n+1 n+1 n+1 n+1 [εk ρk vk φk ]w ⎬ n+1 ⎤ n+1 ⎤ ⎪ ⎡ ⎡ φ φ ∂ ∂ ⎢εk n + 1Λk k ⎥ − ⎢εk n + 1Λk k ⎥ ⎪ ⎢⎣ ⎢⎣ ∂z ⎥⎦ ∂z ⎥⎦ ⎪ e w ⎪ [S(φk )Δz]P + [K (φl n + 1 − φk n + 1)Δz]P ⎭

⎡ n+1 n+1 n+1 Δz ⎤ − εk nρk n φk n) ⎥ ⎢⎣(εk ρk φk Δt ⎦ P + − = +

for k ≠ l

(13)

It is sometimes favorable to linearize the source term according to: S(φk )Δz = (SkC + [Skφφk ]P )Δz

(14)

In general, the linearized discrete equations are presented on the form: aP φP = aW φW + aEφE + SPC =

∑ anbφnb + SPC nb

(15)

where the coefficients aP, aW, and aE depend on the discretization scheme applied. If the differencing scheme produces coefficients that satisfy criterion 16 the resulting matrix of coefficients has the desirable feature of being diagonally dominant.45 This implies that Sφk in eq 14 should always be less than or equal to zero. |aW | + |aE| ⎧≤1 at all nodes ⎨ |aP | ⎩ εmin, otherwise (as)J = 1, (as)J−1 = 0, (as)J+1 = 0, and (bs)J = 0, which sets the solid phase velocity to zero in node J.

■

RESULTS AND DISCUSSION In the present study, a dynamic one-dimensional two-fluid model has been derived. The model has been adopted to simulate the

Figure 8. Void fraction of solid predicted by the two-dimensional twofluid model: 40 s of simulation.

Figure 9. Gas phase velocity predicted by the two-dimensional two-fluid model: 40 s of simulation. Due to symmetry, only half of the domain is computed but copied for complete visualization.

4214



Article

Figure 10. Dry mole fractions at the reactor outlet, i.e. z = H, as a function of time for CO (a), CO2 (b), CH4 (c), and H2 (d). Simulation results of onedimensional model and cross-sectional average of two-dimensional model: SMR and SE-SMR processes.

model results, the following averaging operator is employed on the two-dimensional model results: ⟨ψ ⟩A =

∫A ψ da ∫A da

≈

∑i ψiAi A

(102)

Standard Case Simulation Results. Adopting the reactor operating conditions, bed size, and numerical parameters as specified by Table 7, Figures 4−7 present the comparison between the one- and two-dimensional model simulation results after 40 s of simulation. Considering the SMR process, the one-dimensional model predicts a larger bed expansion than obtained with the twodimensional model, Figure 4a. In Figure 4b, the gas phase velocity compares fairly well comparing the results of the two model approaches in the lower part of the reactor bed, but differences are observed in the upper part due to the relatively larger bed expansion of the one-dimensional model. Further, a significant temperature drop at the reactor entrance is obtained with the one-dimensional model, whereas the simulation results of the two-dimensional model give a uniform temperature, Figure 4c. On the other hand, Figure 5 shows that the species composition compares well between the two model solutions. Adopting the one-dimensional model for simulation of the SESMR process, only small differences are observed in the phase fraction and gas phase velocity predictions relative to that obtained in the simulation of the SMR process, i.e. Figure 6a and

Figure 11. SMR process. Gas phase temperature. Sensitivity to effective heat conductivity, i.e. Cs in eq 105, Table 2, on the simulated temperature profile. One-dimensional model results compared with the standard cross-sectional averaged two-dimensional case.

predict deviations in the variables that differ significantly from reality. The validation of the one-dimensional model are based on the comparison with the corresponding simulation results of a twodimensional model for (i) flow and species composition throughout the reactor bed at 40 s of simulation and (ii) species composition at the reactor outlet as a function of time. In order to allow direct comparison between the one- and two-dimensional 4215



Article

The relative differences in the flow, i.e. εs and vg variables in Figure 4a and b for the SMR process and Figure 6a and b for the SE-SMR process, are related to the CO2 capture, that is, mass of CO2 is removed from the gas phase into the solid phase where CO2 is chemical bounded to the sorbent in accordance to reaction 4. However, due to the effect of mixing the reformer catalyst with the sorbent, this has a larger effect on the simulation results of the two-dimensional model than on the results obtained with the one-dimensional model. The two-dimensional representations in Figures 8 and 9 reveal the significant differences in the flow, i.e. εs and vg, between the SMR and SESMR process. The relatively larger two-dimensional mixing effect observed in the reactor bed for the SMR process compared to the SE-SMR process reveals the reason for the temperature variations in the SE-SMR process, Figure 6c, relative to the homogeneous temperature profile obtained in the SMR process, Figure 4c. Moreover, the one-dimensional model obtains relatively larger temperature drops in the reactor bed compared to the two-dimensional model, Figures 4c and 6c, due to the lack of description of the radial convective flow pattern as obtained by the two-dimensional model, Figures 8 and 9. A dynamic presentation of the species composition of CO, CO2, CH4, and H2 is given in Figure 10. Minor differences are observed in the prediction of the chemical reactor conversion between the one- and two-dimensional models. The differences in the simulation results obtained with the one- and twodimensional models the very first simulation seconds are related to the start-up of the multiphase reactive flow algorithms. Sensitivity to Model ParametersExtended Diffusive and Conductive Fluxes. Due to the differences between the simulation results of the one- and two-dimensional models, sensitivity analysis have been performed on the interphase momentum transfer and extended conductive fluxes. In order to compensate for the radial convective flow pattern, as inherent in the two-dimensional model (Figures 8 and 9), the one-dimensional model may be modified by extended conductive fluxes. The sensitivity to the scaling parameter Cs in eq 105, Table 2, is given in Figure 11 for the SMR process. A relatively high value of Cs is required in order to reduce the large temperature drop at the reactor entrance. For Cs ≈ 105, the temperature profile obtained with the one-dimensional model compares well with the simulation results of the two-dimensional model. The simulation results show negligible sensitivity to the extended diffusive fluxes. This is in agreement with the simulation results presented in Figures 5, 7, and 10, where there are minor differences in the chemical composition predicted by the one- and two-dimensional models. Sensitivity to Model ParametersInterphase Drag. As discussed previously, a larger bed expansion is obtained with the one-dimensional model compared to the bed expansion in the two-dimensional model, Figures 4a and 6a. The bed expansion is correlated to the interphase momentum transfer between the solid and gas phases. Thus if the interphase drag coefficient β in eq 124, Table 5, is, e.g., reduced, the solid particles are less affected by the drag force exerted by the up-flowing gas phase. As a consequent of lowering the interphase momentum transfer between the phases, a reduction in the bed expansion is expected. Figure 12a holds the sensitivity to parameter a in eq 124, Table 5. For a = 35 (a = 150 in the interphase drag correlation suggested by Gidaspow33), the bed expansion obtained with the onedimensional model compares well with the solid fraction profile of the two-dimensional model. On the other hand, a lowering in

Figure 12. SMR process. Sensitivity to interphase momentum transfer on the simulated profiles of (a) εs and (b) vg. Parameter a in eq 124, Table 5. One-dimensional model results compared with the standard cross-sectional averaged two-dimensional case.

Figure 13. SMR process. Sensitivity to interphase momentum transfer on the simulated profile of εs. Parameter b in eq 124, Table 5. Onedimensional model results compared with the standard cross-sectional averaged two-dimensional case.

b versus Figure 4a and b, respectively. On the other hand, relatively large differences are obtained for the εs and vg profiles simulating the SMR and SE-SMR processes with the twodimensional model. Adopting the one-dimensional model, a smaller temperature drop is observed in the temperature prediction within the SE-SMR reactor bed, Figure 6c, relative to that of the temperature prediction of the SMR process, Figure 4c. This observation is due to the heat generated in the reaction between CO2 and the solid sorbent. In contradiction to the twodimensional model result, the temperature drop obtained with the one-dimensional model is concentrated at the reactor entrance and is thus relatively larger. 4216



Article

Figure 14. SE-SMR process: (a) εs, (b) vg, and (c) Tg. Cs = 105 in eq 105, Table 2 and a = 35 in eq 124, Table 5. One-dimensional model results compared with the standard cross-sectional averaged two-dimensional case.

Figure 15. SE-SMR process: (a) H2 and (b) CO2. Cs = 105 in eq 105, Table 2 and a = 35 in eq 124, Table 5. One-dimensional model results compared with the standard cross-sectional averaged two-dimensional case.

the parameter b in eq 124, Table 5, has a minor effect on the bed expansion, Figure 13. Unfortunately, the reduction in the bed expansion achieved by reducing parameter a in the drag coefficient eq 124 results in increased gas phase velocity because

of the denser packing of the solid particles, Figure 12b. For this reason, the one-dimensional model should be considered modified by including the effect of gas bypass in terms of gas bubbles in the dense solid bed. The development of a one4217



Article

dimensional model. For a stronger model validation, the modified one-dimensional model should be evaluated and compared with the two-dimensional model for other reactor operation conditions and reactor bed sizes.

dimensional model with the incorporation of the gas bubble effect should result in a lowering in the gas phase velocity, and as consequence, the driving force of the interphase momentum term should be reduced. Figure 14 holds the reactor bed profiles of εg, vg, and Tg and Figure 15 gives the species fractions of H2 and CO2. The onedimensional model is modified by an extended conductive flux corresponding to Cs = 105 in eq 105, Table 2, and modified drag correlation corresponding to a = 35 in eq 124, Table 5. Good agreement is obtained between the one- and two-dimensional models for the species fractions. Moreover, fair results can be obtained for εg and Tg with the one-dimensional model. On the other hand, the well-fitted bed expansion to the two-dimensional result by modification of the interphase drag coefficient is at the cost of the increased differences in the gas phase velocity prediction of the one- and two-dimensional models. Limitations of the Current Numerical Solution Algorithm. There are several challenges related to the numerical solution of the dynamic one-dimensional two-fluid model. As with the multidimensional two-fluid models, the numerical solution algorithm of the current model is sensitive to the strong pressure−velocity coupling. Further, for multiphase reactive gas−solid flows extra effort is placed on the numerical solution algorithm due to the gas phase density variations caused by the chemical reactions. The modeling of a bubbling fluidized bed reactor also cause problems associated with the gradients coming from the transition zone between the dense solid bed and the free board zone. Moreover, the finite volume method adopted in this study is a low order method with a low convergence rate. Computational Costs. The 40 s of simulations of the SMR and SE-SMR processes by the two-dimensional two-fluid model implemented in FORTRAN takes 3−4 days on a supercomputer; whereas, the corresponding one-dimensional simulation results are obtained within 5−6 h on a standard personal computer using the programming language MATLAB. However, the programming language MATLAB was used for the simulation of the onedimensional model due to the convenient debugging and graphical properties. Translating the MATLAB language into, e.g., FORTRAN will reduce the computational time which will be of particular interests when the model is extended from a stand alone fluidized bed to a circulating fluidized bed design with interconnected solid fluxes transferred between the reformer and regenerator units. Further Work. The Kunii−Levenspiel type of models are empirical to a large extent. Thus, fair simulation results may be obtained fitting the parameters of the Kunii−Levenspiel type of models to the simulation results obtained with the twodimensional model. Because the absence of phasic momentum coupling terms in the Kunii−Levenspiel type of models, the problems associated with the bed expansion and lack of gas bypassing are a model limitation of the one-dimensional model. The Kunii−Levenspiel model and the one-dimensional two-fluid models should be evaluated and compared modeling the SESMR process in the circulating fluidized bed design. In such a study the importance of the inherent dynamic solid flux in the two-fluid model, compared to the prescribed solid flux in the Kunii−Levenspiel type of models, will be elucidated. A well-defined interphase drag coefficient β is important in order to describe the gas−solid flow within the reactor bed. Thus, appropriate drag formulation should be developed for onedimensional models. Moreover, modification of the model including the gas bypassing phenomena by bubbles may be a solution to the overpredicted gas phase velocity of the one-

■

CONCLUSION A dynamic one-dimensional two-fluid model has been derived and applied to simulated the SMR and SE-SMR processes operated in bubbling fluidized bed reactors. The simulation results have been compared with the corresponding results of a two-dimensional model. The simulation results of the one-dimensional model are in good agreement with the two-dimensional model considering the chemical conversion of the reactor, as is also utilized by the Kunii−Levenspiel type of models. Moreover, with extended conductive fluxes, fair temperature profiles can be predicted with the one-dimensional model. On the other hand, the flow pattern, i.e. the phasic fractions and gas phase velocity, is associated with the largest uncertainty in current model. However, the internal flow details do not have significant influence on the chemical process performance. Nevertheless, if the present two-fluid model can be extended by including the effect of gas bypassing, improved description of the flow will presumably be obtained. However, external solid flow and heat transfer may still be important. Thus, the current one-dimensional model is considered to have good potential for further model development in order to study interconnected fluidized bed reactors with a dynamic solid flux transferred between the reactor units.

■

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected] (J.S.); hugo. [email protected] (H.A.J.). Notes

The authors declare no competing financial interest.

■ ■

ACKNOWLEDGMENTS The fellowship (J.S.) financed by the Research Council of Norway (NFR) is gratefully appreciated. NOMENCLATURE

Latin Letters

⟨φ⟩A = area averaging quantity A = problem coefficient matrix of the discretization a = coefficient of the discretization a′ = coefficient of the discretization Cp = heat capacity, J/(kg K) ) = coefficient used in PEA; heat balances b = coefficient of the discretization; source term b′ = coefficient of the discretization; source term dk = velocity correction coefficient of phase k D = mass dispersion coefficient Dir = binary diffusion coefficient, m2/s + = coefficient used in PEA; momentum balances Dmg,i = Wilke diffusion coefficient, m2/s Dk = submatrix of matrix A holding the coupling terms d = diameter, m - = coefficient used in PEA; generalized f f = fanning friction factor F = convective flux, kg/(m2 s) g = gas gz = gravity, m/s2

4218



Article

G = modulus of elasticity, kg/(s2 m) ΔH = heat of reaction, J/kmol H = height of reactor bed, m Hmf = height of solid in the reactor bed at minimum fluidization, m hgs = interfacial heat transfer coefficient, kg/(s3 K) I = generic flux k = generalized phase index, k = g,s k = heat conductivity, J/(m s K) K = generalized proportional coefficient in interphase flux term l = generalized phase index, l ≠ k M = molecular mass, kg/kmol max(a,b) = maximum of a,b 4k = interphase momentum transfer, kg/(m2 s2) pg = gas phase pressure, Pa pcorrection = gas phase pressure correction, Pa g p*g = guessed gas phase pressure, Pa Qi = interfacial heat transfer, J/(m3s) R = gas constant, J/(kmol K) Rrxi = reaction rate of reaction rxi, kmol/(m3 s) Rk,i = formation rate of species i in phase k, kmol/(m3 s) Re = Reynolds number, dimensionless t = time, s Δt = time step, s T = temperature, K s = solid S = generalized source term SC = source term independent of the dependent variable Sφ = source term coefficient dependent on the dependent variable sign(a) = returns the sign of a vk = velocity of phase k, m/s vcorrection = velocity correction of phase k, m/s k vcorrected = corrected phase k velocity due to pressure, m/s k z = axial reactor dimension, m Δz = axial resolution, m

k = phase, k = g,s l = phase, l ≠ k mf = minimum fluidization nb = neighbor node rxi = reaction index s = solid t = tube, reactor bed P = central node p = solid particle wall = reactor bed wall w = west interfacial node W = west node Superscript

Cap = CO2 capture reaction eff = effective i = interface in = reactor entrance m = molecular n = time level out = at the reactor outlet SMR = reaction of the steam methane reforming stat = static t = turbulent ν = iteration Abbreviation

■

PEA = partial elimination algorithm 1D = one-dimensional 2D = two-dimensional SMR = steam methane reforming SE-SMR = sorption-enhanced steam methane reforming

REFERENCES

(1) Balasubramanian, B.; Lopez Ortiz, A.; Kaytakoglu, S.; Harrison, D. P. Hydrogen from methane in a single-step process. Chem. Eng. Sci. 1999, 54, 3543−3552. (2) Rostrup-Nielsen, J. R. Catalytic steam reforming. In Catalysis science and technology; Anderson, J. R., Boudant, M., Eds.; Springer: Berlin, 1984; Vol. 4, pp 1−118. (3) Williams, R. Hydrogen production. U.S. Patent 1,938,202, 1933. (4) Abanandes, J. C.; Rubin, E. S.; Anthony, E. J. Sorbent cost and performance in CO2 capture systems. Ind. Eng. Chem. Res. 2004, 43, 3462−3466. (5) Li, Z.-S.; Cai, N.-s. Modeling of multiple cycles for sorptionenhanced steam methane reforming and sorbent regeneration in fixed bed reactor. Energy Fuels 2007, 21, 2909−2918. (6) Ochoa-Fernández, E. CO2 acceptors for sorption-enhanced steam methane reforming. Ph.D. thesis, Norwegian University of Science and Technology (NTNU), 2007. (7) Rusten, H. K.; Ochoa-Ferández, E.; Chen, D.; Jakobsen, H. A. Numerical investigation of sorption enhanced steam methane reforming using Li2ZrO3 as CO2-acceptor. Ind. Eng. Chem. Res. 2007, 46, 4435− 4443. (8) Gupta, H.; Fan, L. S. Carbonation−calcination cycle using high reactivity calcium oxide for carbon dioxide separation from flue. Ind. Eng. Chem. Res. 2002, 41, 4035−4042. (9) Sun, P.; Grace, J. R.; Lim, C. J.; Anthony, E. J. Determination of intrinsic rate constants of the CaO−CO2 reaction. Chem. Eng. Sci. 2008, 63, 47−56. (10) Alvarez, D.; Abanades, J. C. Determination of the critical product layer thickness in the reaction of CaO with CO2. Ind. Eng. Chem. Res. 2005, 44, 5608−5615. (11) Mess, D.; Sarofim, A. F.; Longwell, J. P. Product layer diffusion during the reaction of calcium oxide with carbon dioxide. Energy Fuels 1999, 13, 999−1005.

Greek Letters

Γk = interphase mass transfer flux, kg/(m3 s) εk = cross-sectional area fraction of phase k ε′g = area fraction of gas phase that is not in terms of bubbles εbubbles = area fraction of gas phase in bubbles g ρ = density, kg/m3 ρ*, ρ**, ρ*** = density estimates by fractional step τ = shear stress, kg/(s2 m) ψ = generalized quantity ω = species mass fraction φ = generalized property Λ = generalized diffusion coefficient εiter = relative iteration error εres = residual error γ = Defined by eq 39, kg/(K s3 m) μ = viscosity, kg/(m s) β = interphase drag coefficient Subscript

CO2 = carbon dioxide e = east interfacial node E = east node g = gas i = species type j = grid index, velocity node J = grid index, scalar node 4219



Article

(12) Stendardo, S.; Foscolo, P. U. Carbon dioxide capture with dolomite: A model for gas−solid reaction within the grains of a particulate sorbent. Chem. Eng. Sci. 2009, 64, 2343−2352. (13) Solsvik, J.; Jakobsen, H. A. A numerical study of a two property catalyst/sorbent pellet design for the sorption-enhanced steam-methane reforming process: Modeling complexity and parameter sensitivity study. Chem. Eng. J. 2011, 178, 407−422. (14) Abanades, J. C. The maximum capture efficiency of CO2 using a carbonation/calcination cycle of CaO/CaCO3. Chem. Eng. J. 2002, 90, 303−306. (15) Alvarez, D.; Abanades, J. C. Pore-size and shape effects on the recarbonation performance of calcium oxide submitted to repeated calcination/recarbonation cycles. Energy Fuels 2005, 19, 270−278. (16) Grasa, G. S.; Abanades, J. C. CO2 capture capacity of CaO in long series of carbonation/calcination cycles. Ind. Eng. Chem. Res. 2006, 45, 8846−8851. (17) Fernandez, J. R.; Abanades, J. C.; Murillo, R. Modeling of sorption enhanced steam methane reforming in an adiabatic fixed bed. Chem. Eng. Sci. 2012, 84, 1−11. (18) Koumpouras, G. C.; Alpay, E.; Stepanek, F. Mathematical modelling of low-temperature hydrogen production with in situ CO2 capture. Chem. Eng. Sci. 2007, 62, 2833−2841. (19) Rout, K. R. A study of the sorption-enhanced steam methane reforming process. Ph.D. thesis, Norwegian University of Science and Technology (NTNU), 2012. (20) Wang, Y.; Chao, Z.; Jakobsen, H. A. Numerical study of hydrogen production by the sorption−enhanced steam methane reforming process with online CO2 capture as operated in fluidized bed reactors. Clean Technol. Environ. Pol. 2011, 13, 559−565. (21) Wang, Y.; Chao, Z.; Chen, D.; Jakobsen, H. A. SE-SMR process performance in CFB reactors: Simulation of the CO2 adsorption/ desorption process with CaO based sorbents. Int. J. Greenhouse Gas Control 2011, 5, 489−497. (22) Wang, Y.; Chao, Z.; Jakobsen, H. A. 3D simulation of bubbling fluidized bed reactors for sorption enhanced steam methane reforming processes. J. Nat. Gas Sci. Eng. 2010, 2, 105−113. (23) Lindborg, H.; Jakobsen, H. A. Sorption enhanced steam methane reforming process performance and bubbling fluidized bed reactor design analysis by use of a two−fluid model. Ind. Eng. Chem. Res. 2009, 48, 1332−1342. (24) Kunii, D.; Levenspiel, O. Fluidization Engineering; ButterworthHeinemann: Boston, 1991. (25) Kunii, D.; Levenspiel, O. Fluidized reactor models. Ind. Eng. Chem. Resour. 1990, 29, 1226−1234. (26) Kunii, D.; Levenspiel, O. Circulating fluidized bed reactors. Chem. Eng. Sci. 1997, 52, 2471−2482. (27) Kunii, D.; Levenspiel, O. The K−L reactor model for circulating fluidized beds. Chem. Eng. Sci. 2000, 55, 4563−4570. (28) Shimizu, T.; Satoh, M.; Fujikawa, F.; Tonsho, M.; Inagaki, M. Simultaneous reduction of SO2, NOx, and N2O emissions from a twostage bubbling fluidized bed combustor. Energy Fuels 2000, 14, 862− 868. (29) Kunii, D.; Levenspiel, O. Bubbling bed model for kinetic processes in fluidized beds. Ind. Eng. Chem. Process Design Dev. 1968, 7, 481−492. (30) Romano, M. C. Modeling the carbonator of a Ca-looping process for CO2 capture from power plant flue gas. Chem. Eng. Sci. 2012, 69, 257−269. (31) Jakobsen, H. A. Chemical reactor modeling: Multiphase reactive flows; Springer-Verlag: Berlin Heidelberg, 2008. (32) Lun, C. K. K.; Savage, S. B.; Jeffrey, D. J. Kinetic theories for granular flow: inelastic particles in Couette flow and slightly inelastic particles in general flow field. J. Fluid Mech. 1984, 140, 223−256. (33) Gidaspow, D. Multiphase flow and fluidization. Continuum and kinetic theory descriptions; Academic Press, New York, 1994. (34) Lindborg, H.; Lysberg, M.; Jakobsen, H. A. Practical validation of the two-fluid model applied to dense gas−solid flows in fluidized beds. Chem. Eng. Sci. 2007, 62, 5854−5869.

(35) Wang, Y.; Chao, Z.; Jakobsen, H. A. A sensitivity study of the twofluid model closure parameters (β, e) determining the main gas−solid flow pattern characteristics. Ind. Eng. Chem. Res. 2010, 49, 3433−3441. (36) Satrio, J. A.; Shanks, B. H.; Wheelock, T. D. Development of a novel combined catalyst and sorbent for hydrocarbon reforming. Ind. Eng. Chem. Res. 2005, 44, 3901−3911. (37) Satrio, J. A.; Shanks, B. H.; Wheelock, T. D. A combined catalyst and sorbent for enhanced hydrogen production from coal or biomass. Energy Fuels 2007, 21, 322−326. (38) Albrecht, K. O.; Satrio, J. A.; Shanks, B. H.; Wheelock, T. D. Application of a combined catalyst and sorbent for steam reforming of methane. Ind. Eng. Chem. Res. 2010, 49, 4091−4098. (39) Chao, Z.; Wang, Y.; Jakobsen, J. P.; Fernandino, M.; Jakobsen, H. A. Derivation and validation of a binary multi-fluid Eulerian model for fluidized beds. Chem. Eng. Sci. 2011, 66, 3605−3616. (40) Carlo, A. D.; Bocci, E.; Zuccari, F.; Dell’Era, A. Numerical investigation of sorption enhanced steam methane reforming process using computational fluid dynamics Eulerian−Eulerian code. Ing. Eng. Chem. Res. 2010, 49, 1561−1576. (41) Chao, Z.; Wang, Y.; Jakobsen, J. P.; Fernandino , M.; Jakobsen, H. A. Numerical investigation of the sorption enhanced steam methane reforming in a fluidized bed reactor. Energy Procedia 2012, 26, 15−21. (42) Rout, K. R.; Solsvik, J.; Nayak, A. K.; Jakobsen, H. A. A numerical study of multicomponent mass diffusion and convection in porous pellets for the sorption-enhanced steam methane reforming and desorption processes. Chem. Eng. Sci. 2011, 66, 4111−4126. (43) Xu, J.; Froment, G. F. Methane steam reforming, methanation and water−gas shift: I. Intrinsic kinetics. AIChE J. 1989, 35, 88−96. (44) Patankar, S. V. Numerical heat transfer and fluid flow; Hemisphere Publishing Corporation: London, 1980. (45) Versteeg, H. K.; Malalasekera,W. An introduction to computational fluid dynamics: The finite volume method, 2nd ed.; Pearson Prentice Hall: Upper Saddle River, NJ, 2007. (46) Ferziger, J. H.; Períc, M. Computational methods for fluid dynamics, 3rd ed.; Springer: Berlin, 2002. (47) Spalding, D. B. Turbulent buoyant convection; Hemisphere: Washington, 1977; pp 569−586, ICHMT seminar 1976. (48) Spalding, D. B. In Recent advances in numerical methods in fluids; Taylor, C., Morgan, K., Eds.; Pineridge Press: Swansea, 1980; Vol. 1, pp 139−167. (49) Spalding, D. B. Mathematical methods in nuclear-reactor thermal hydraulics. In Proceedings of ANS meeting on nuclear reactor thermal hydraulics; Lahey, R. T., Ed.; Saratoga NY, Oct 8, 1980; pp 1979−2023. (50) Spalding, D. B. In Computational techniques in heat transfer; Lewis, R., Morgan, K., Johnson, J. A., Smith, W. R., Eds.; Pineridge Press, 1985; pp 1−44. (51) Spalding, D. B. New development and computed results; Imperial College of Science and Technology, London, 1981. (52) Carver, M. B. Numerical computation of phase separation in two fluid flow. J. Fluids Eng. 1984, 106, 147−153. (53) Lindborg, H. Modeling and simulations of reactive two-phase flows in fluidized beds. Ph.D. thesis, The Norwegian University of Technology and Science (NTNU), 2008. (54) Jakobsen, H. A. On the modeling and simulation of bubble column reactors using a two-fluid model. Ph.D. thesis, The Norwegian Institute of Technology (NTH), Trondheim, Norway, 1993. (55) Chao, Z.; Wang, Y.; Jakobsen, J. P.; Fernandino, M.; Jakobsen, H. A. Investigation of the particle−particle drag in a dense binary fludized bed. Powder Technol. 2012, 224, 311−322. (56) Wilke, C. R. Diffusional properties of multicomponent gases. Chem. Eng. Prog. 1950, 46, 95−104. (57) Reid, R. C.; Prausnitz, J. M.; Poling, B. E. The properties of gases and liquids; McGraw-Hill: New York, 1987. (58) Wen, C. Y.; Yu, Y. H. Mechanics of fluidization. Chem. Eng. Prog. Symp. Ser. 1966, 62, 100−111.

4220


Sorption-Enhanced Steam

Recommend Documents