Parametric study of hydrogen production via sorption enhanced steam methane reforming in a circulating fluidized bed riser

Chemical Engineering Science 192 (2018) 1041–1057

Contents lists available at ScienceDirect

Chemical Engineering Science journal homepage: www.elsevier.com/locate/ces

Parametric study of hydrogen production via sorption enhanced steam methane reforming in a circulating fluidized bed riser Kiattikhoon Phuakpunk a, Benjapon Chalermsinsuwan b,c,d,⇑, Sompong Putivisutisak d,e, Suttichai Assabumrungrat a,d a

Center of Excellence in Catalysis and Catalytic Reaction Engineering, Department of Chemical Engineering, Faculty of Engineering, Chulalongkorn University, Bangkok 10330, Thailand Fuels Research Center, Department of Chemical Technology, Faculty of Science, Chulalongkorn University, Bangkok 10330, Thailand Center of Excellence on Petrochemical and Materials Technology, Chulalongkorn University, Bangkok 10330, Thailand d Advanced Computational Fluid Dynamics Research Unit, Chulalongkorn University, Bangkok 10330, Thailand e Department of Mechanical Engineering, Faculty of Engineering, Chulalongkorn University, Bangkok 10330, Thailand b c

h i g h l i g h t s
Kinetics of steam methane reforming and carbonation were incorporated in CFD model. k


2 factorial design was used to analyze parameters.
Hydrogen purity could reach equilibrium as well as that conducted in bubbling bed.
Gas velocity, riser diameter and solid flux had major effects on system designing.

a r t i c l e

i n f o

Article history: Received 31 January 2018 Received in revised form 13 August 2018 Accepted 18 August 2018 Available online 22 August 2018 Keywords: Sorption enhanced steam methane reforming Computational fluid dynamics Circulating fluidized bed Multiphase flow models Riser

a b s t r a c t Computational fluid dynamics was applied for sorption enhanced steam methane reforming (SESMR) operating in a circulating fluidized bed (CFB) riser. The solid mixtures consisted of Ni-based catalyst and CaO sorbent. The aim of study was to design a proper pilot-scale CFB riser which produced hydrogen (H2) with both high purity and high flux. The design parameters and the reaction parameters were examined with 2k full factorial design. The significances of each parameter were analyzed by analysis of variance. Using the optimum result, the highest H2 purity reached 98.58% in dry basis accompanied with the highest H2 flux of 0.301 kg/m2 s. The hydrodynamics of this optimum case showed that SESMR was nearly completed since 5.0 m height because axial and radial distributions of solid were well developed without excessive segregation between catalyst and sorbent. Thus, the H2 purity and the H2 flux approached fully developed within the riser height. Ó 2018 Elsevier Ltd. All rights reserved.

DH 298 ¼ þ206:2 kJ=mol

1. Introduction

bSMR :

CH4 þ H2 O
CO þ 3H2

Conventional industrial-level hydrogen production has used methane (CH4), in natural gas or in tail gas from refinery process, as raw material via steam methane reforming (SMR) process. The conventional processes consist of three main sections i.e. reforming, shifting and gas separation (Harrison, 2008). In reforming and shifting section, maximum CH4 is converted to mainly CO2 and H2 (with little of CO) via three reversible reactions including steam methane reforming (bSMR and gSMR) reactions and water-gas shift (WGS) reaction as shown below.

gSMR :

CH4 þ 2H2 O
CO2 þ 4H2

⇑ Corresponding author at: Fuels Research Center, Department of Chemical Technology, Faculty of Science, Chulalongkorn University, Bangkok 10330, Thailand. E-mail address: [email protected] (B. Chalermsinsuwan). https://doi.org/10.1016/j.ces.2018.08.042 0009-2509/Ó 2018 Elsevier Ltd. All rights reserved.

ð1Þ

DH 298 ¼ þ165:0 kJ=mol ð2Þ

WGS :

CO þ H2 O
CO2 þ H2

DH 298 ¼ 41:2 kJ=mol

ð3Þ

By thermodynamic equilibrium, a high temperature above 750 °C is sufficient to maximize conversion of CH4 in a reforming furnace. The effluent gas of the furnace which still containing CO about 8–10% in dry basis is then fed to WGS reactors including high temperature shift (HTS) reactor, operating at 300–400 °C, and low temperature shift (LTS) reactor, operating at 200–300 °C, in series. The outlet gas from the second shift reactor approximately consists of 76% H2, 17% CO2, 4% unreformed CH4 and 3% CO in dry basis

1042

K. Phuakpunk et al. / Chemical Engineering Science 192 (2018) 1041–1057

Nomenclature CD Cfr;ls Cat/Sb ds ess  els f drag Gs * g g0;ls g0;ss H Hg Hj Hfj hsg  hgs I id

* J qk

Ki Kk Kls  Ksl Ksg  Kgs kcarb kg ki kHs Mj _ pq m Nus n Pr p pCO2 pCO2 ;eq pg

Drag coefficient, [–] Friction coefficient between solid phases, [–] Catalyst to sorbent ratio, [kgkg1] Particle diameter of solid phase, [m] Restitution coefficient for solid-solid collisions, [–] Drag function Solid flux, [kgm2s1] Gravity force, [ms2] Radial distribution coefficient of mutual solid phases, [–] Radial distribution coefficient of single solid phase, [–] Height of the riser, [m] Specific enthalpy of gas phase, [m2s2, Jkg1] Specific enthalpy of species j in the reaction, [m2s2, Jkg1] Specific heat of formation of species j in the reaction, [m2s2, Jkg1] Gas-solid interphase heat exchange coefficient, [m2s2, Jkg1] Unit tensor, [–] Diameter of the riser, [m] Mass flux of species k into phase q, [kgm2s1] Equilibrium constants of reaction i Adsorption equilibrium constants of species k Solid-solid interphase momentum exchange coefficient, [kgm3s1] Gas-solid interphase momentum exchange coefficient, [kgm3s1] Rate constants of carbonation Thermal conductivity of gas phase, [Wm1K1] Rate constants of reaction i Diffusion coefficient, [m2s1] Molecular weight of species j in the reaction, [kgkmol1] Mass transfer from phase p to phase q, [kgm3s1] Nusselt number of solid phase, [–] Degree of partial pressure, [–] Prandtl number of gas phase, [–] Static pressure, [Pa] Partial pressure of CO2, [Pa] Equilibrium pressures of CO2, [Pa] Static pressure of gas phase, [Pa]

(Harrison, 2008). Thereafter, H2 is separated through separation units such as pressure swing adsorption (PSA) or amine scrubbing technology (Barelli et al., 2008). These conventional processes use many units, utilities and resources. A new concept of sorption enhanced steam methane reforming (SESMR) has been raised because of the advantages taken from that CO2 that adsorbed by sorbents in the same reforming reactor. For the first advantage, equilibrium of reforming reactions (1) and (2) is shifted forwardly, so H2 is more produced and CH4 can be almost completely converted. Another advantage is that effluent gas from the reformer has higher H2 purity reaching of 99% in dry basis (Harrison, 2008; Barelli et al., 2008; Koumpouras et al., 2007; Cotton et al., 2013), thus separation units are unnecessary in the processes. However, limitation of SESMR is the discontinuous performance when the sorbent is almost full of CO2 captured. The sorbent has to release CO2 before reprocessing SESMR. There are many sorbents suitable for CO2 capture divided into natural sorbents such as limestone, dolomite, huntite and hydrotalcite, and synthetic sorbents such as Li4SiO4, Li2ZrO3 and Na2ZrO3

pk ps Q sg * qg

 R Rqk

Res rcarb ri Sh;g Sqk Sm;q S0 S/C Tin U * v ls * vq * v rj * v sg

XCaO Yqk

cj eq es;max Hs

kq ks

lq ls ls;col ls;fr ls;kin qq sq ss

Partial pressures of species k, [Pa] Solid pressure, [Pa] Intensity of heat exchange between solid phase and gas phase, [Wm3] Heat flux of gas phase, [Wm2] Rate of a heterogeneous reaction, [kmolm3s1] Net rate of species k produced by homogeneous reactions inside phase q, [kgm3s1] Particle Reynolds number of solid phase, [–] 1 Rate of carbonation, [kmolkg1 ] sorbs 1 Rates of reaction i, [kmolkg1 s ] cat Heat source of gas phase, [Wm3] Rate of creation of species k by addition from dispersed phase and other sources in phase q, [kgm3s1] Mass source of phase q, [kgm3s1] Initial specific surface area of CaO, [m2kgsorb 1] Steam to carbon, [molmol1] Temperature of inlets, [°C] Gas inlet velocity, [ms1] Interphase velocity from phase l (solid or gas) to solid phase, [ms1] Velocity of phase q, [ms1] Velocity of a reactant j which involved in the reaction, [ms1] Interphase velocity from solid phase to gas phase, [ms1] Conversion of CaO, [–] Mass fraction of species k in phase q, [–] Stoichiometric coefficient of species j in the reaction, [–] Volume fraction of phase q, [–] Maximum packing of solid phase, [–] Granular temperature, [m2s2, Jkg1] Bulk viscosity of phase q, [Pas] Solid bulk viscosity, [Pas] Shear viscosity of phase q, [Pas] Solid shear viscosity, [Pas] Collisional viscosity of solid phase, [Pas] Friction viscosity of solid phase, [Pas] Kinetic viscosity of solid phase, [Pas] Physical density of phase q, [kgm3] Stress tensor of phase q, [Pa] Particulate relaxation time in solid phase, [s]

(Harrison, 2008; Barelli et al., 2008). Among all listed sorbents, CaO sorbents, i.e. limestone and dolomite, are the cheapest and have the highest CO2 capacity, moderate adsorption rate but the lowest stability (Harrison, 2008; Barelli et al., 2008). Hydrotalcite-like materials (HTCls) have the highest adsorption rate, good stability but very low CO2 capacity (Harrison, 2008; Barelli et al., 2008; Koumpouras et al., 2007). The synthetic sorbents have high capacity, good stability and low to moderate adsorption rate but are the most expensive (Ochoa-Fernandez et al., 2005; Harrison, 2008; Barelli et al., 2008; Koumpouras et al., 2007). Thus, for large-scale hydrogen production via SESMR, natural CaO sorbents such as dolomite and limestone are preferable due to their profitable costs. Even though limestone has more CO2 capacity than dolomite due to higher CaO content, but dolomite contains more MgO that makes dolomite more stable to cyclic usage (Harrison, 2008; Comas et al., 2004; Aceves Olivas et al., 2014). Thus, dolomite is a suitable sorbent for large-scale processes of SESR which has to involve regeneration of the used sorbent.

K. Phuakpunk et al. / Chemical Engineering Science 192 (2018) 1041–1057

A forward reaction of CaO adsorbing CO2 is called carbonation and its backward reaction is called decarbonation as the following reaction.

CaO þ CO2
CaCO3

DH 298 ¼ 178 kJ=mol

ð4Þ

Fluidized bed reactors are widely applied in processes containing gas and solids. Advantages of fluidized bed reactors over fixed bed reactors are lower pressure drop and more gas-solid (reactantcatalyst/sorbents) contacting throughout beds that give better mass and heat transfers. Bubbling bed reactor has been proven to give SESMR performance as good as the fixed bed reactors (Chao et al., 2012; Wang et al., 2010; Johnsen et al., 2006a, 2006b). Circulating fluidized bed reactors (CFBR) have been utilized in commercial processes such as combustion boiler, fluidized bed catalytic cracking (FCC), Fischer-Tropsch synthesis, etc. For large hydrogen production via SESMR, CFBR has been proposed to use because CFBR gives higher productivity and lower pressure drop than a bubbling bed reactor. In CFBR, feed and produced gas can be increased by using higher gas inlet velocity of the reformer. Another advantage of CFBR is that the reformer is unnecessary to be switched to regeneration mode because sorbents can be transported out of the reformer and regenerated simultaneously in another reactor (Kunii and Levenspiel, 1991). However, higher gas velocity, gas-solid contacting is theoretically non-uniform among dilute zone and dense zone. From these reasons, CFBR has unclear performance which is dependent on flow regimes (Arstad et al., 2012; Wang et al., 2011; Rodríguez et al., 2011), so the CFBR system design has been interesting for SESMR application with high efficiency. In the last decade, computational fluid dynamics (CFD) has attracted more interest and has been utilized in problems about fluid phenomena especially in system geometry design. CFD could give high accuracy and detailed results that could reduce time, resources and costs of experiments especially in large scales. Problem solving can be either in 2D or 3D models, the results from 3D models are more close to realistic but consume much more computational demands and time than from 2D models. However, in many cases especially in cylindrical geometry like the riser, 2D models could demonstrate sufficient accuracy details (Wu et al., 2013, 2014; Wang et al., 2014; Di Carlo et al., 2010; Samruamphianskun et al., 2012; Chalermsinsuwan et al., 2014). Johnsen et al. (2006a) extended the SESMR experiment, in a single bubbling bed reformer, by scaling-up to pilot scale and included a regenerator in dual bubbling fluidized bed reactors as well as Sánchez et al. (2012). Johnsen et al. (2006a) simulated the system with constant circulating solids but Sánchez et al. (2012) simulated the system with solid return between the reactors, simultaneously. The results were sufficient for more practical design and prediction of realistic phenomena inside CFBR system by non-commercial CFD programs. However, the dual bubbling fluidized bed reactors systems from both Johnsen et al. (2006a) and Sánchez et al. (2012) were quite complicated for practical operation. Thus, the reformer in this study became a riser, like in a typical CFBR, with higher gas velocity for blowing solids out in practice. Moreover, a commercial CFD program like ANSYSÒ FluentÒ could help users to perform simulations of fluidized bed reactors combined with complex kinetics of SESMR especially when scaling up or considering modified configuration of the reformer like in a work of Herce et al. (2017) as well as in this study. Herce et al. (2017) used ANSYSÒ FluentÒ to upsize a bubbling bed reactor from good-performance laboratory-scale to the large scale. They studied the effect of temperature on SESMR performance. They found that a modified drag model could predict more accurate results for the large scale bubbling bed simulation. However, the major aims of this study were to investigate the feasibility of performing SESMR in a riser reactor for higher hydrogen produc-

1043

tion and to study parametric sensitivity for further experimental development in the future. In this study, ANSYSÒ FluentÒ 15.0.7 with 2D models was performed to design a pilot-scale riser suitable for the best performance of hydrogen production via SESMR by investigating effects of design and reaction parameters including riser diameter, inlet temperature, catalyst to sorbent ratio, solid flux and gas velocity. The performance of the system was represented by H2 flux and H2 purity of the effluent gas and both were response variables for 2k full factorial design analysis. The optimum case was predicted and additionally simulated to analyze hydrodynamics in the system. 2. Methodology 2.1. System description According to Kunii and Levenspiel (1991) and Prajongkan (2011), a basic design of CFBR system is shown in Fig. 1 including a reformer, a cyclone and a regenerator. The reformer could be a riser due to very fast kinetics of reforming and carbonation reactions whereas the regenerator needed sufficient residence time for slower kinetics of decarbonation to nearly complete removing undesired CO2. Thus, the regenerator should be a bubbling bed reactor or a semi-batch reactor. The key part of this CFBR system was the riser, and thus design and reaction parameters of the riser would be studied to determine feasibility of SESMR performing in CFBR. The simulated results from chosen values of the parameters had to satisfy constraints/goals as follows: 1. There is no solid accumulation in the riser. 2. The effluent gas achieved high H2 purity at 99% in dry basis (theoretical equilibrium), while H2 production (H2 flux) was as high as possible. 3. CaO conversion (XCaO) in solid outlet should be lower than 28% which was the stable maximum capacity of dolomites after steady re-cycling, according to Johnsen et al. (2006a). Together with the results of Johnsen et al. (2006b) and Sánchez et al. (2012), a new preliminary design of the reforming riser was typically calculated as in Kunii and Levenspiel (1991). The riser has sufficient 7 m height and around 0.1 m diameter with inlet of fed gas at the bottom, inlet of solid at a side channel and outlet of the mixture at the top. The channel of inlet solid has 0.05 m diameter with its center is at 0.075 m high above the bottom. The solid included two individual phase, i.e. catalyst phase and sorbent phase. The SMR reactions would be in the catalyst phase while the carbonation would be in the sorbent phase. Properties of each solid are presented in Table 1. The catalyst used in simulation represented Ni-based catalyst and was assumed to have uniform size and density. The fed sorbent was the fresh dolomite whose properties were assumed uniform as well. The fresh dolomite contained 60% weight of CaO and the other of MgO. 2.2. Computation description and models In the last decade, there have been two main approaches for gas-solid flow modeling, i.e. Euler-Euler model and EulerLagrange model (Ranade, 2002; Yeoh and Tu, 2010). In the EulerEuler model, also called granular flow model (GFM), an Eulerian framework is considered for all phases. But in the Euler-Lagrange model, also called discrete particle models (DPM), the Eulerian framework is considered for the continuum phase (gas phase) and a Lagrangian framework is considered for all the dispersed phases (solid phases). The Euler-Euler model has been more suitable for simulating complex and large amount of solid particles like

1044

K. Phuakpunk et al. / Chemical Engineering Science 192 (2018) 1041–1057

Fig. 1. The system configuration and parameters of the CFBR riser.

The rate expressions for the reactions (1)–(3) are provided in Eqs. (5)–(7), respectively.

Table 1 The conditions and properties in models.

Catalyst density Calcined dolomite density Mean catalyst particle size Mean dolomite particle size MgO content in dolomite Packing limit of catalyst and sorbent phase Restitution coefficient of all phase interactions

2200 1540 200 250 40 0.60 0.90

kg/m3 kg/m3 lm lm wt% – –

fluidized bed (Chao et al., 2012; Wang et al., 2010, 2011, 2014; Di Carlo et al., 2010; Sánchez et al., 2012; Prajongkan, 2011; Sánchez et al., 2012, 2013; Solsvik et al., 2012, 2014; Wang et al., 2010). The concept of Euler-Euler model is that summation of volume fraction of each phase in each control volume is always unity and the solid phases are considered as continuum as fluid phase. The kinetic theory of granular flow (KTGF) was then additionally applied to solid phases. The KTFG has extended from kinetic theory of gas by adding kinetic energy oscillation owing to inelastic collisions and fluctuating motions of the particles (Gidaspow, 1994). A representative of the kinetic energy oscillation comes in terms of a granular temperature (Hs ) and can be evaluated from kinetic fluctuation energy conservation. Hydrodynamic models which mainly consisted of conservation equations of mass, momentum, energy and species together with many constitutive equations were performed for flow phenomena in domains. The heterogeneous kinetic equations of SMR and carbonation, which had been used in previous researches (Chao et al., 2012; Wang et al., 2010, 2011; Sánchez and Jakobsen, 2012; Sánchez et al., 2013; Solsvik et al., 2012, 2014; Wang et al., 2010), were simultaneously employed in volumetric form. All concerning models are presented below. 2.2.1. Kinetics of SESMR Xu and Froment (1989) investigated the intrinsic kinetics which accounted for the resistance to diffusion on a Ni/MgAl2O4 catalyst.

!

pCH4 pH2 O 

rgSMR ¼

kgSMR p3:5 H2

pCH4 p2H2 O 

rWGS ¼

  pH pCO2 kWGS 1  pCO pH2 O  2 pH2 KWGS DEN2

rbSMR

Conditions and properties

p3H2 pCO

kbSMR ¼ 2:5 pH2

KbSMR p4H2 pCO2 KgSMR



1

ð5Þ

DEN2

! 

1

ð6Þ

DEN2

DEN ¼ 1 þ KCO PCO þ KH2 PH2 þ KCH4 PCH4 þ

ð7Þ KH2 O PH2 O PH2

ð8Þ

Rate constants of reaction i are expressed by Xu and Froment (1989) as follows:    240100 1 1 1 kbSMR ¼ 9:708  104 exp  ; ½kmol  Pa0:5  kgcat  s1  R T 648 ð9Þ kgSMR ¼ 1:156  104 exp



  243900 1 1 1  ; ½kmol  Pa0:5  kgcat  s1  R T 648 ð10Þ

kWGS ¼ 1:2597  106 exp



  67130 1 1 1  ; ½kmol  kgcat  Pa1  s1  R T 648 ð11Þ

According to Xiu et al. (2002), equilibrium constants of reaction i are as follows: KbSMR ¼

1013252 expð0:2513Z 4  0:3665Z 3  0:58101Z 2 þ 27:1337Z  3:277Þ

;

½Pa2 

ð12Þ

KgSMR ¼ KbSMR KWGS ;

½Pa2 

ð13Þ

1045

K. Phuakpunk et al. / Chemical Engineering Science 192 (2018) 1041–1057

KWGS ¼ expð0:29353Z 3 þ 0:63508Z 2 þ 4:1778Z þ 0:31688Þ;

½ ð14Þ

* @q þ r  ðq v Þ ¼ Sm @t

ð25Þ

For each phase q:

1000  1; Z¼ T

½

ð15Þ

Xu and Froment (1989) also expressed adsorption equilibrium constants of species k as follows:

KCH4

   38280 1 1  ; ¼ 1:791  10 exp R T 823

KH 2 O

   88680 1 1 ¼ 0:4152 exp  ; R T 823

6

1

½Pa 

½

ð16Þ

ð17Þ

   70650 1 1  ; ¼ 4:091  10 exp R T 648

½Pa 

ð18Þ

   82900 1 1 KH2 ¼ 2:960  107 exp  ; R T 648

½Pa1 

ð19Þ

KCO

4

1

Sun et al. (2008) determined the rate constants of carbonation of two calcined natural CaO sorbents i.e. limestone and dolomite. The kinetics of carbonation expressed below is applied for dolomite.

rcarb ¼ kcarb ðpCO2  pCO2 ;eq Þn S0 ð1  XCaO Þ

ð20Þ

Equilibrium pressure of CO2 (pCO2 ;eq ) is dependent on temperature ranges as follows:

n X * @ _ pq  m _ qp Þ þ Sm;q ðeq qq Þ þ r  ðeq qq v q Þ ¼ ðm @t p¼1

b) Momentum conservation The conservation equation of momentum can be written as * * * * @ ðq v Þ þ r  ðq v v Þ ¼ rpþr  s þ q g @t

pCO2 ;eq

* * * * @ ðeg qg v g Þ þ r  ðeg qg v g v g Þ ¼ eg rpþr  sg þ eg qg g @t n X * * * * _ sg v sg  m _ gs v gs Þ þ ðKsg ðv s  v g Þ þ m s¼1

ð28Þ And the momentum conservation for solid phase (s) is * * * * @ ðes qs v s Þ þ r  ðes qs v s v s Þ ¼ es rp  rps þ r  ss þ es qs g @t n X * * * * _ ls v ls  m _ sl v sl Þ þ ðKls ðv l  v s Þ þ m l¼1

ð29Þ ð21Þ

- For T 6 1; 173:15 K according to Johnsen et al. (2006);

  20474 pCO2 ;eq ¼ 4:1918  1012 exp T

ð22Þ

Rate constants of carbonation (kcarb ) and degree of partial pressure (n) are dependent on pCO2 ranges as follows: - For ðpCO2  pCO2 ;eq Þ > 10; 000 Pa;

kcarb

  20400 ; ¼ 1:04  106 exp RT

½kmol  m2  s1 

ð27Þ

But unlike mass conservation, interaction between solid phase and solid phase is different from interaction between fluid phase and solid phase. Thus, for gas phase (g) the momentum conservation can be

- For T > 1,173.15 K according to Abanades et al. (2004);

  19130 ¼ 1:216  1012 exp T

ð26Þ

ð23Þ

where solid pressure (ps), stress tensor of phase q (sq ) and interphase momentum exchange coefficient (Ksg ; Kls ) are described in the constitutive equations.

c) Energy conservation Analogous to velocities in momentum conservation, enthalpies take places in terms of mass flow in conservation of energy. For gas phase the equation of energy conservation can be written as a below.

@pg * * * @ ðeg qg Hg Þ þ r  ðeg qg v g Hg Þ ¼ eg þ sg : r v g  r  q g @t @t n X _ sg hsg  m _ gs hgs Þ þ Sh;g þ ðQ sg þ m s¼1

ð30Þ

n¼0 - For 0 < ðpCO2  pCO2 ;eq Þ 6 10; 000 Pa;

  20400 kcarb ¼ 1:04  1010 exp ; RT

½kmol  m2  Pa1  s1 

ð24Þ

n¼1 - For ðpCO2  pCO2 ;eq Þ 6 0 Pa, there is no CO2 captured i.e. decarbonation occurs instead.

where interphase heat exchange coefficient between solid phase and gas phase (hsg ; hgs ) is described in constitutive equations. For a solid phase, the energy conservation is derived with KTGF and expressed as the kinetic fluctuation energy conservation as follows:

  * * 3 @ ðes qs Hs Þ þ r  ðes qs v s Hs Þ ¼ ðps I þ ss Þ : r v s 2 @t

þ r  ðkHs rHs Þ  cHs þ /ls ð31Þ *

2.2.2. Governing equations The governing equations are not only conservative equations of mass, momentum and energy but also conservation of chemical species and their source terms due to heterogeneous reactions.

where ðps I þ ss Þ : r v s term is generation of energy by the solid stress tensor, r  ðkHs rHs Þ term is diffusion of energy. cHs is collisional dissipation of energy, defined by Lun et al. (ANSYS Inc., 2013) as

a) Mass conservation The conservative equation of mass or called ‘‘continuity equation” can be expressed as

cHs ¼

12ð1  e2ss Þg0;ss pffiffiffiffi qs e2s Hs3=2 ds p

ð32Þ

1046

K. Phuakpunk et al. / Chemical Engineering Science 192 (2018) 1041–1057

/ls is kinetic energy exchange between phase l (could be either solid or gas) and solid phase, defined as

/ls ¼ 3Kls Hs

ð33Þ

where granular temperature (Hs ) and radial distribution coefficient (g0;ss ) are described in constitutive equations. d) Chemical species conservation General equation of chemical species conservation for a species k in a single phase is in the following form. * * @ ðqYk Þ þ r  ðq v Yk Þ ¼ r  J k þ Rk þ Sk @t

e) Radial distribution coefficient For one solid phase the radial distribution coefficient (g0;ss ) is modified by Lun et al. (ANSYS Inc., 2013) as a below.

"



g0;ss ¼ 1 

es

1=3 #1 ð42Þ

es;max

For mutual solid phases between phase s and another solid phase l, the mutual radial distribution coefficient (g0;ls ) is given as

g0;ls ¼

ds g0;ll þ dl g0;ss ds þ dl

ð43Þ

ð34Þ

In multiphase flow, the conservation of the species k in phase q is *q *q @ q q q ðe q Yk Þ þ r  ðeq qq v Yqk Þ ¼ r  ðeq J k Þ þ eq Rqk þ eq Sqk @t n X _ pk qk  m _ qk pk Þ þ ðm

ð35Þ

f) Granular temperature from KTGF The granular temperature (Hs ) is the representative of energy in KTGF as described previously and can be evaluated from the kinetic fluctuation energy conservation as Eq. (31). Algebraic formulation is able to simplify the conservation equation by neglecting convection and diffusion terms as the following equation. * 3 @ ðes qs Hs Þ ¼ ðps I þ ss Þ : r v s  cHs þ /ls 2 @t

p¼1

ð44Þ

2.2.3. Constitutive equations a) Stress tensor In multiphase, the stress tensor for phase q is *

*



2 3



*

sq ¼ eq lq ðr v q þ r v Tq Þ þ eq kq  lq r  v q I

g) Gas-solid momentum exchange coefficient The gas-solid momentum exchange coefficient (Ksg ) is defined

ð36Þ

Ksg  Kgs ¼

For solid phase, shear viscosity and bulk viscosity are described next. b) Solid shear viscosity In solid phase the shear viscosity involves viscosity due to collision, kinetics and friction (optional). Thus the solid shear viscosity (ls ) is given as

ls ¼ ls;col þ ls:kin þ ls;fr

as

ð37Þ

ls;col

ss ¼

qs d2s 18lg

and the kinetic viscosity is applied from Gidaspow et al. (ANSYS Inc., 2013) as

ls;kin

 2 10qs ds Hs p 4 ¼ 1 þ es g0;ss ð1 þ ess Þ es 5 96es ð1 þ ess Þg0;ss

ks ¼

 1=2 4 2 H es qs ds g0;ss ð1 þ ess Þ s 3 p

ps ¼ es qs Hs þ 2e qs Hs g0;ss ð1 þ ess Þ

*

24

eg Res

- For

½1 þ 0:15ðeg Res Þ0:687 

ð48Þ

eg 6 0:8;

Ksg ¼ 150

ð47Þ

*

*

es ð1  eg Þlg qg es jv s  v g j þ 1:75 2 ds eg ds

ð49Þ

ð40Þ h) Solid-solid momentum exchange coefficient The solid-solid momentum exchange coefficient (Kls , l is another solid phase) is expressed as

d) Solid pressure The solid pressure (ps ) as a function of the granular temperature was given by Lun et al. (ANSYS Inc., 2013). The expression consisting of kinetic energy and particle collision terms is 2 s

*

ð39Þ

c) Solid bulk viscosity The solid bulk viscosity (ks ) expression by Lun et al. (ANSYS Inc., 2013) is

eg > 0:8;

3 es eg qg jv s  v g j Ksg ¼ CD eg  2:65 ds 4 CD ¼

pffiffiffiffiffiffiffiffiffiffi

ð46Þ

Due to drag function (f drag ) is dependent on drag coefficient (CD ) and one of the proper drag functions is Gidaspow’s drag model (ANSYS Inc., 2013), the gas-solid momentum exchange coefficients from Gidaspow’s drag model are as follows: - For

ð38Þ

ð45Þ

where particulate relaxation time in solid phase (ss ) is given as

where the collisional viscosity is

 1=2 4 Hs ¼ es qs ds g0;ss ð1 þ ess Þ es 5 p

es qs f drag ss

ð41Þ

Kls  Ksl ¼

  2 2 3ð1 þ els Þ p2 þ Cfr;ls p8 es qs el ql ðdl þ ds Þ g0;ls 2pðql dl þ qs ds Þ 3

3

*

*

jv l  v s j

ð50Þ

K. Phuakpunk et al. / Chemical Engineering Science 192 (2018) 1041–1057

i) Gas-solid heat exchange coefficient The interphase heat exchange coefficient between solid phase and gas phase (hsg ) is expressed as

hsg  hgs ¼

kg Nus ds

ð51Þ

where Nusselt number of solid phase (Nus) has a correlation according to Gunn’s model (ANSYS Inc., 2013) as 1=3 Nus ¼ ð7  10eg þ 5e2g Þð1 þ 0:7Re0:2 Þ þ ð1:33  2:4eg s Pr 1=3 þ 1:2e2g ÞRe0:7 s Pr

ð52Þ

The computational domain of the each riser configuration as shown in Fig. 2 was drawn by ANSYSÒ DesignModelerTM and its mesh was created by ANSYSÒ MeshingTM. The uniform domains were chosen with 4 different cell sizes of mesh refinement. The following simulations were calculated using ANSYSÒ FluentÒ in a series of ANSYSÒ products (version 15.0.7). The SESMR operation was simulated for 20 s using the time step of 1  103 s. At initial time, there was only inert N2 in the domain. Then both fed gas and solid had entered through their inlet boundaries in normal direction. The calculated constant velocities and volume fraction were the inputs of each phase but the inlet granular temperature was set for only solid phases. The fed gas was fixed to consist of CH4 and H2O in a ratio of 1:4 by mole (S/C = 4) which had been investigated suitable and sufficient for SMR (Joensen and Rostrup-Nielsen, 2002). All feeds came in with the same temperature while the wall was adiabatic.

1047

The wall surface was set in no-slip condition. The product mixture discharges the riser to atmosphere. All parameters and properties of the system to be fixed or studied are shown in Table 2.

2.3. Parametric study From Fig. 1, there are a large number of parameters including the design parameters (sizes of reformer, gas velocity and solid flux) and the reaction parameters (reaction temperature, gas composition and solid composition). All of these parameters could mutually affect the performance of H2. Because one-factor-at-atime (OFAT) investigating would take lots of simulation cases to find the best case. A statistic method like 2k full factorial design, which chose only 2 expected levels in each of k parameters, could make only 2k different independent cases to investigate. This method was also able to analyze significance of each single parameter and their interactions via analysis of variance (ANOVA) method. This method has linear assumption between the two experimental points because it is a screening method. Only five parameters of the system including the gas inlet velocity, the solid flux, the diameter of riser, the catalyst to sorbent ratio and the temperature of inlets were explored with 25 factorial design. Two levels of these parameters and other system properties are shown in Table 2. The response variables i.e. H2 flux and H2 purity represented the performance of H2 production. Finally, the optimized case of SESMR in the riser would be determined.

Fig. 2. The computational domains of the riser with different cell sizes.

1048

K. Phuakpunk et al. / Chemical Engineering Science 192 (2018) 1041–1057

Table 2 The parameters and system properties chosen in the 25 full factorial design. Parameters

Value

Design parameters Gas inlet velocity (U) Solid flux (Gs) Diameter of the riser (id) Height of the riser (H) Reaction parameters Catalyst to sorbent ratio (Cat/Sb) Steam to carbon (S/C) Temperature of inlets (Tin) CaO conversion of inlet sorbent (XCaO,in) Outlet pressure Inlet granular temperature Wall type Shear condition

4 and 6 50 and 200 0.05 and 0.2 7

m/s kg/m2 s m m

0.16 and 2.54 4 575 and 665 0 101,325 1  105 Adiabatic No slip

kg/kg mol/mol °C % Pa m2/s2

3. Results and discussion 3.1. Validations

tions at different gas velocities are shown in Fig. 3. At the lowest gas velocity (0.320 m/s), several small bubbles occurred and dispersed. When gas velocity was increased, some small bubbles collapsed into bigger bubbles and the bed was lifted higher. With much higher gas velocity (up to 0.892 m/s), the bubbles were much larger and more expanded the bed height. The time-averaged bed heights of the bed from these simulations, experiments of Lin et al. (1985) and simulations of Sánchez et al. (2012) are summarized in Table 3. The results indicated that all of bed heights were very close. In SESMR validation, UDFs of SESMR kinetics were compiled to validate with experimental results from the bubbling bed reformer by Johnsen et al. (2006b). This validation was separated into two cases, at first only SMR kinetics was compiled then additional carbonation kinetics was employed for SESMR in the other case. The comparison of simulation effluent gas compositions to experimental ones are shown in Table 4. Most of them had good agreement with slight deviation. Different Ni content, different structure and different properties might affect the CO and CO2 compositions which were sensitive with WGS reaction.

There were two validations performed in this study i.e. cold flow validation and SESMR validation. In cold flow validation, chemical reactions were not involved in simulations. This validation was performed to test all chosen hydrodynamic models. The cold flow validation compared the bed heights of bubbling bed to reference experiments of Lin et al. (1985) and simulation of Sánchez et al. (2012). The instantaneous contour of solid volume fraction and velocity vectors of solid at 10 s resulted from simula-

Because instantaneous flow of each particle and amount of each substance had fluctuated within fluidized bed reactors, thus, timedependent simulation had to be performed rather than steadystate simulation. A quasi-steady state value from selected timeaveraged range was used to represent the result. The riser with

(a) U = 0.320 m/s

(b) U = 0.458 m/s

(c) U = 0.641 m/s

3.2. Time average and mesh refinement

(d) U = 0.892 m/s

Fig. 3. The instantaneous solid volume fraction (left) and solid velocity (middle) of the bubbling bed reactor at 10 s relating to the experimental solid mean velocity (right) of Lin et al. (1985) with various gas velocities.

1049

K. Phuakpunk et al. / Chemical Engineering Science 192 (2018) 1041–1057 Table 3 The comparison of time-averaged bed heights of cold flow validation with the experimental results of Lin et al. (1985) and the simulation results of Sánchez et al. (2012).

U = 0.320 m/s U = 0.458 m/s U = 0.641 m/s U = 0.892 m/s

Lin et al.

Sánchez et al.

This simulation

0.12 0.15 0.21 0.23

0.145 0.16 0.185 0.225

0.130 0.155 0.182 0.234

0.3

H2 flux [kg/m2.s]

Bed height [m]

0.35

0.25 0.2 0.15 0.1 0.05

3.3. The 25 factorial design analysis for H2 production

0 0.0

2.0

4.0

6.0

8.0

10.0

12.0

14.0

16.0

18.0

20.0

Time [s] Fig. 4. The H2 flux out as a function of time in case of id = 0.1 m, Tin = 600 °C, Cat/ Sb = 0.16, U = 6 m/s and Gs = 200 kg/m2 s.

7

6

5

4

Height [m]

diameter of 0.1 m operating with inlet temperature of 600 °C, catalyst to sorbent ratio of 0.16, gas velocity of 6 m/s and 200 kg/m2 s solid flux, H2 flux out of the riser at every 0.1 s was plotted as an example case in Fig. 4. The results showed that the fluctuation of H2 flux seemed stable after approximately 5 s as well as in other cases. Thus, a time-averaged range of 10–20 s would be used to represent further simulation results in this study. Next, to investigate the mesh refinement, the same example cases were performed and time-averaged axial profiles of H2 flux with different sizes of cell (DxDy) were displayed as shown in Fig. 5. The results showed that the 10 mm20 mm size was insufficiently fine while the 5 mm10 mm size took more unnecessary calculating time. Although the 5 mm20 mm size and the 10 mm10 mm size had equivalent the numbers of cells but the 5 mm20 mm size was chosen for further simulation because in the riser, radial profile had more effects than the axial profile.

3

2

First of all, an example simulation case would be visually analyzed. The instantaneous contour plots of volume fraction of solid phases (catalyst and sorbent) and gas phase at 20 s in case with the riser had diameter of 0.2 m operated with inlet temperature of 575 °C, catalyst to sorbent ratio of 2.54, gas velocity of 6 m/s and 200 kg/m2 s solid flux (run No. 20) are displayed in Fig. 6. The gradient shades of colors from blue to red represent low volume fraction to high volume fraction of each phases. The contour plots show that both the catalyst and the sorbent phases were denser near the wall all along the height of the riser while the volume fraction of gas phase was dense in the center line of the riser. These volume fraction contour plots confirmed that the fast fluidization occurred in the riser. The similarity of volume fraction of the catalyst phase and the sorbent phase indicated there was no segregation between catalyst and sorbent. This made good mixing of the catalyst and the sorbent, thus SESMR was performed very well in this CFB riser. Fig. 7 displays the instantaneous contour plot of mole fraction of H2 in the gas phase at 20 s in the same case of run No. 20. The contour plot shows the development of H2 produced along the height of the riser. In this shown case, H2 was completely in equilibrium before exiting the riser. Area-averaged fractions of H2 and H2O near the outlet were 0.6387 and 0.3504, respectively, and the little remaining fraction was the other gases. This fraction of H2 equaled to H2 purity of 98.37% in dry basis.

10mm x 20mm 10mm x 10mm 1

5mm x 20mm 5mm x 10mm

0 0

0.05

0.1

0.15

0.2

0.25

0.3

0.35

Time-averaged axial H2 flux [kg/m2.s] Fig. 5. The time-averaged axial profiles of H2 flux with different cell sizes in case of id = 0.1 m, Tin = 600 °C, Cat/Sb = 0.16, U = 6 m/s and Gs = 200 kg/m2 s.

Considering temperature contour plots of all phases in run No. 20 as shown in Fig. 8, temperature of each phase had the same pattern that temperature changed slightly from initial 848 K (575 °C) to minimum 833 K in left zone, to maximum 854 K in middle zone and rarely changed in right zone. Because total SMR reactions are highly endothermic (DH°298 = +165.0 kJ/mol of CH4), while carbonation is highly exothermic (DH°298 = 178 kJ/mol of CH4), therefore SESMR is slightly exothermic (DH°298 = 13 kJ/mol of CH4). The upward flows with almost constant temperature on the right side indicated that SESMR performed well (SMR and carbonation performed in balance) since solid entrance. The upward flows with temperature decrease on the left side showed that SMR on the catalyst occurred more than carbonation on the sorbent (both SMR and carbonation performed but SMR was some higher). This confirmed by

Table 4 The comparison of SMR and SESMR validation with experimental results of Johnsen et al. (2006b).

SMR SESMR

Johnsen et al. This simulation Johnsen et al. This simulation

H2 [% dry]

CH4 [% dry]

CO [% dry]

CO2 [% dry]

73.4 74.3 98 97.7

6 6.4 1 1.6

8 2.7 0.5 0.4

12 16.6 0.5 0.3

1050

K. Phuakpunk et al. / Chemical Engineering Science 192 (2018) 1041–1057

(a) Catalyst

(b) Sorbent

(c) Gas

Fig. 6. The instantaneous volume fraction of each phase at 20 s of run No. 20.

the H2 fraction which developed higher on left side of the bottom zone as in Fig. 7. Considering the increase of temperature accompanied with further development until full of the H2 fraction in the middle zone, this could described that CO2 from the left zone was adsorbed and made SMR shifted forwardly (both SMR and carbonation performed but carbonation was some higher). Next to analyze 25 factorial design results from all 32 runs, the area average near the outlet with the time average of 10–20 s was used to report values of each H2 flux, H2 purity and others, e.g. CaO conversion, as shown in Table 5. The H2 flux and the H2 purity, which both represented H2 production performance, were response variables of 32 runs of the factorial design. As previously described constraint, Johnsen et al. (2006a) specified that the limitation of using dolomite as the sorbent in CFBR was losing CaO capacity in long term re-usage until the CaO conversion could not be over 28% steadily. In every run, none of CaO conversion reached 28%, thus the limitation of CaO capacity was not needed to be concerned and long-term circulating of the dolomite was feasible in ranges of this studied system. The other results showed that there was no run which get the maximum H2 flux together with the maximum H2 purity, thus both needed further statistical analysis like the ANOVA for determining the best practice and sensitivity analyses. The riser diameter (id), the inlet temperature (Tin), the catalyst to sorbent ratio (Cat/Sb), the solid flux (Gs), the gas velocity (U), the H2 flux and the H2 purity were coded as A, B, C, D, E, R1 and R2, respectively, in the ANOVA. In the ANOVA, any main effect or interaction which significantly affected the H2 flux or the H2 purity had to have P-value less than

0.05 (Montgomery, 2012). The ANOVA results of both the H2 flux and the H2 purity are shown in Tables 6 and 7, respectively. The main effects and their interactions in significant order were descending sorted by obvious to the p-values. Furthermore, regression models for prediction the H2 flux and the H2 purity were determined as Eqs. (53) and (54), respectively.

H2 flux ¼ 0:222599 þ 0:019A þ 0:015B þ 0:014C þ 0:017D þ 0:030E þ 0:00970AE þ 0:00967DE

ð53Þ

H2 purity ¼ 93:97 þ 2:55A þ 1:04C þ 2:43D  2:71E  1:36AB  2:00AD þ 1:48AE  1:21BD þ 1:41DE  1:19ADE ð54Þ where A, B, C, D and E in the regression models were coded variables which transformed low to high levels of considered parameters into 1 to +1 range. The P-values from both the tables and the coefficients in the both regression models indicated that most top-three significant main effects/interactions on both the H2 flux and the H2 purity were the same, i.e. the gas velocity (E), the riser diameter (A) and the solid flux (D) in descending order. The next descending significant main effects/interactions were different between on the H2 flux and on the H2 purity. The inlet temperature (B) significantly affected only the H2 flux but its interactions (AB and BD) occurred significant to the H2 purity. On overall H2 production performance, the design parameters including the gas velocity (E), the riser

1051

K. Phuakpunk et al. / Chemical Engineering Science 192 (2018) 1041–1057

(a) Catalyst

(b) Sorbent

(c) Gas

Fig. 8. The instantaneous temperature of each phase at 20 s of run No. 20.

Fig. 7. The instantaneous mole fraction of H2 (wet basis) in gas phase at 20 s of run No. 20.

diameter (A) and the solid flux (D), in descending order, had more effect than the reaction parameters including the inlet temperature (B) and the catalyst to sorbent ratio (C). Next for the sensitivity analyses, the main effects and interactions on the H2 flux and the H2 purity were plotted as shown in Figs. 9 and 10, respectively. The slopes showed either positive or negative effects on the response variable and the steepness of slope could indicate the significant order like the P-value and the regression coefficients. In case that the H2 flux was the response variable in Fig. 9a, all of the main effects were positive as well as their coefficients in the regression model. Among these main effects, the slope of the gas velocity (E) was obviously the steepest, following by the riser diameter (A), the solid flux (D), the inlet temperature (B) and the catalyst to sorbent ratio (C), in descending order. The positive effect of the gas velocity (E) described that increasing of the gas velocity made feed increased despite the residence time must be less. But for the rapid reaction like SESMR, the reaction had sufficient time to produce more H2. In case of the riser diameter (A), the larger diameter get more H2 flux because the bigger area of gas inlet made the feed increased as well, even though the dilute region of solid might be wider. Comparison of the solid flux (D) indicated that the low flux (50 kg/m2 s) was insufficient contact between gas (reactant) and solid (catalyst/sorbent). Lastly, both of the inlet temperature (B) and the catalyst to sorbent ratio (C) prefer the high level (665 °C and 2.54 kg/kg, respectively) comparing to the low levels.

From the interaction effects on the H2 flux as plotted in Fig. 9b, when the gas velocity (E) was operated with the low riser diameter (A) or the low solid flux (D), the slopes had the same direction as the main effect (E) but were slightly less steep. On the other hand, when the gas velocity (E) was operated with the high riser diameter (A+) or the high solid flux (D+), the slopes had the same direction as the main effect (E) but were slightly steeper. Because all of the main effects (E, A and D) were positive, thus their interactions (AE and DE) were more positive to the H2 flux. The high gas velocity interacted with the high riser diameter (A+E+) get higher H2 flux because they both made more feed mutually. However, the little changes of the steepness indicated that the interactions (AE and DE) were less effective than the main effect (E) and corresponding to their P-values in Table 6 and regression coefficients in Eq. (53). In Fig. 10a, the riser diameter (A), the solid flux (D) and the catalyst to sorbent ratio (C) had positive effects on the H2 purity as well as H2 flux. This meant that at the high levels of these parameters (id = 0.2 m, Gs = 200 kg/m2 s and Cat/Sb = 2.54 kg/kg), contact of gas–solid (reactant-catalyst/sorbent) was better than at the low levels. Only the gas velocity (E) had negative effect on the H2 purity which was opposed to the H2 flux. This is because the higher gas velocity (6 m/s) made more feed but the residence time was less and insufficient to reach SESMR equilibrium. When the gas velocity (E) was operated with the low riser diameter (A) or the low solid flux (D) in Fig. 10b, the slopes of the H2 purity had the same negative direction as the main effect (E) but were slightly steeper. Whereas when the gas velocity (E) was operated with the high riser diameter (A+) or the high solid flux (D+), the slopes had the same negative direction as the main effect (E) but were slightly less steep. This meant that even though

1052

K. Phuakpunk et al. / Chemical Engineering Science 192 (2018) 1041–1057

Table 5 The H2 flux, the H2 purity and the CaO conversion (XCaO) from simulations with the 25 factorial design. Factor: Run

A id [m]

B Tin [°C]

C Cat/Sb [kg/kg]

D Gs [kg/m2s]

E U [m/s]

R1 H2 flux [kg/m2s]

R2 H2 purity [% dry]

XCaO [%]

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32

0.05 0.05 0.05 0.05 0.05 0.05 0.05 0.05 0.05 0.05 0.05 0.05 0.05 0.05 0.05 0.05 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2

575 575 575 575 575 575 575 575 665 665 665 665 665 665 665 665 575 575 575 575 575 575 575 575 665 665 665 665 665 665 665 665

2.54 2.54 2.54 2.54 0.16 0.16 0.16 0.16 2.54 2.54 2.54 2.54 0.16 0.16 0.16 0.16 2.54 2.54 2.54 2.54 0.16 0.16 0.16 0.16 2.54 2.54 2.54 2.54 0.16 0.16 0.16 0.16

50 50 200 200 50 50 200 200 50 50 200 200 50 50 200 200 50 50 200 200 50 50 200 200 50 50 200 200 50 50 200 200

4 6 4 6 4 6 4 6 4 6 4 6 4 6 4 6 4 6 4 6 4 6 4 6 4 6 4 6 4 6 4 6

0.172422 0.189089 0.197854 0.267603 0.13552 0.109459 0.178654 0.21108 0.201309 0.258555 0.206029 0.303812 0.176828 0.177717 0.201382 0.272869 0.200761 0.27316 0.204256 0.297392 0.190798 0.216396 0.203383 0.279256 0.206551 0.305924 0.203183 0.297275 0.202109 0.276525 0.205234 0.300787

93.90 79.71 98.36 94.86 89.73 70.37 96.43 90.63 96.80 88.17 97.67 96.56 94.66 82.56 97.34 95.03 98.60 95.32 99.17 98.37 97.47 91.04 99.10 97.48 97.72 96.93 96.08 91.89 96.36 95.26 97.51 95.98

12.70 10.04 4.01 5.28 3.37 2.55 1.25 1.41 14.94 14.92 3.89 5.85 4.38 4.04 1.33 1.77 4.02 5.32 1.07 1.50 1.27 1.45 0.34 0.43 3.92 5.50 0.92 1.08 1.19 1.71 0.31 0.40

Table 6 The results of the ANOVA of the H2 flux. Source

Sum of squares

Degree of freedom (DF)

Mean square

F-value

P-value

E (U) A (id) D (Gs) B (Tin) C (Cat/Sb) AE DE Residual Cor total

0.02824 0.011356 0.009009 0.006874 0.006249 0.003008 0.002995 0.013187 0.080918

1 1 1 1 1 1 1 24 31

0.02824 0.011356 0.009009 0.006874 0.006249 0.003008 0.002995 0.000549

51.39508 20.66646 16.39574 12.51014 11.37263 5.474424 5.45038

<0.0001 0.000132 0.000465 0.001681 0.002523 0.027946 0.028259

Table 7 The results of the ANOVA of the H2 purity. Source

Sum of squares

Degree of freedom (DF)

Mean square

F-value

P-value

E (U) A (id) D (Gs) AD AE DE AB BD ADE C (Cat/Sb) Residual Cor Total

235.049 207.6093 189.5516 128.5148 69.79685 63.32517 59.46191 46.79862 45.29247 34.3478 125.8546 1205.602

1 1 1 1 1 1 1 1 1 1 21 31

235.049 207.6093 189.5516 128.5148 69.79685 63.32517 59.46191 46.79862 45.29247 34.3478 5.993074

39.2201 34.64154 31.62844 21.44388 11.64625 10.56639 9.921771 7.808784 7.557468 5.731249

<0.0001 <0.0001 <0.0001 0.000144 0.002619 0.003826 0.004833 0.010865 0.012023 0.026079

the higher riser diameter or the higher solid flux made better contact of gas-solid enhancing the H2 purity, the less residence time with the higher gas velocity still had more effect to lower the H2 purity. This was consistent to ANOVA Table 7 which the gas veloc-

ity (E) was much more significant (higher P-value) than their interaction (AE and DE). Similarly in Fig. 10c, when the riser diameter (A) was operated with the inlet temperature (B) or the solid flux (D), the slopes were

1053

K. Phuakpunk et al. / Chemical Engineering Science 192 (2018) 1041–1057

Interacons

H2 flux [kg/m2.s]

Main Effects 0.29

A-Riser id.

0.27

B-Inlet temp.

0.29

(a)

C-Cat./Sorb.

0.25

0.27 0.25

D-Gs 0.23

A-

(b)

A+ DD+

0.23

E-U

0.21

0.21

0.19

0.19 0.17

0.17 -1

-1

1

1

A, B, C, D, E

E

Fig. 9. The main effects and the interactions on the H2 flux.

Main Effects

Interacons

98%

98%

(a)

H2 purity [% dry]

96% 94%

(b)

96% 94%

92%

92%

A-

A-Riser id. 90%

90%

A+

C-Cat./Sorb. 88%

D-Gs

88%

86%

E-U

86%

DD+

-1

1

-1

A, C, D, E

E

Interacons

Interacons

98%

98%

(c)

96%

H2 purity [% dry]

1

94%

(d)

96% 94%

92%

B-

92%

90%

B+

90%

BD-

88%

B+

88%

D+ 86%

86% -1

1

-1

A

1

D

Fig. 10. The main effects and the interactions on the H2 purity.

the same in positive direction as the main effect (A) and slightly steeper at the low level (B and D) but slightly less steep at the high level (B+ and D+). As well as in Fig. 10d, when the solid flux (D) interacted with the inlet temperature (B), the slopes showed the same in positive direction as the main effect (D) with being slightly steeper at the low inlet temperature (B) but being slightly less steep at the high inlet temperature (B+). Lastly, from all of the interaction plots in Fig. 10, the little changes of the steepness of the interactions (AE, DE, AB, AD and BD) indicated that these interactions were less effective than the main effects (A, D and E) consistent with their P-values in Table 7 and regression coefficients in Eq. (54).

Table 8 shows the optimized results of the regression models of the H2 flux and the H2 purity. Either the H2 purity or the H2 flux was optimized prior to the another, the predicted H2 purity was maximum at 99.17% in dry basis and the H2 flux reached 0.309 kg/m2 s in the riser with 0.2 m diameter, the inlet temperature of 581 °C, the catalyst to sorbent ratio of 2.54 kg/kg, the solid flux of 200 kg/m2 s and the gas velocity of 6 m/s. To confirm the results from the prediction, another case with those values of parameters was simulated as run 33. The H2 purity and the H2 flux from the simulation were 98.58% in dry basis and 0.301 kg/m2 s, respectively, which both agreed very well with the results from the prediction.

1054

K. Phuakpunk et al. / Chemical Engineering Science 192 (2018) 1041–1057

Table 8 The optimum H2 flux and H2 purity predicted from the regression models and a simulation. Factor

A id [m]

B Tin [°C]

C Cat/Sb [kg/kg]

D Gs [kg/m2 s]

E U [m/s]

R1 H2 flux [kg/m2 s]

R2 H2 purity [% dry]

Optimizing: R2 prior to R1 Optimizing: R1 prior to R2 Run 33

0.200 0.200 0.2

581.48 581.40 581

2.540 2.540 2.54

199.98 199.82 200

6.00 6.00 6

0.308708 0.308523 0.300574

99.17 99.17 98.58

3.4. Hydrodynamics in the riser In this part, the catalyst and sorbent distributions along axial and radial direction were expected key hydrodynamics which explained how SESMR performed well in this riser system. Inside the riser of the optimum case (run 33), the volume fraction of both catalyst and sorbent phase, volumetric catalyst to sorbent ratio were compared to the H2 flux and the H2 purity. In Fig. 11, the time-averaged of the H2 flux in axial direction was plotted from 0.1 m height which was the exact height above the solid inlet channel. The H2 flux accumulated rapidly in the lower height and reached steady in the upper height near the outlet (7 m). Fig. 12 shows the time-averaged radial distributions of the H2 purity at different heights. Radial distance at 0.0 m was the position of left wall and at 0.2 m was the position of right wall where solid inlet was on this right side at 0.05–0.1 m height. At above solid inlet (0.1 m height), the H2 purity suddenly approached equilibrium at the solid inlet (0.2 m distance). The H2 purity decreased along to left direction far from the solid inlet but was higher at the left wall. At higher height, the H2 purity profile was rapidly higher even in the middle distance. Until over 5.0 m height, the H2 purity profiles approached the equilibrium and were quite steady and uniform along the radial distance. From both Figs. 11 and 12, these profiles could explain that SESMR was extremely rapid. It approached equilibrium since solid (catalyst/sorbent) had started to contact the gas (reactant) at the inlet. In this system of run 33, SESMR could be close to complete since 5.0 m height thus the 7 m height of the riser was sufficient to design.

Because SESMR was a very rapid reaction, solid distribution and mixing of the catalyst and the sorbent were expected as the keys to make SESMR developed in the system as previous discussion. Fig. 13 shows time- and area-averaged axial profiles of volume fraction of the solid phases. The amount of catalyst, sorbent and total solid were most dense at height over the solid inlet (volume fraction of total solid was about 0.25). All the volume fractions of each solid were gradually decreased along axial direction until the volume fraction of total solid was in 0.06–0.20 which was the range of fast fluidized bed regime and close to be in 0.01– 0.06 range of pneumatic transport regime (Kunii and Levenspiel, 1997). Due to dense zone in the lower height thus the axial profile of H2 flux in Fig. 11 increased obviously in this lower zone. In Fig. 14, the time- and area-averaged axial profile of volumetric catalyst to sorbent ratio was plotted comparing with the ratio where the solid entranced (1.778 vol./vol.). The results showed that in the lower zone, the ratio was lower but still more than 1.7 vol./vol. then back to the inlet ratio in the zone above 5.0 m height. This meant there was no serious segregation between the catalyst and the sorbent, and thus SESMR could be well performed all along axial direction. The ratio decreased because the catalyst had higher density (2200 kg/m3) than the sorbent (approximately 1540 kg/m3) so the lighter sorbent was lifted easier.

100%

90%

80%

Time-averaged radial H2 purity [%]

7

6

Height [m]

5

4

3

2

70%

60% 0.1m 1.0m

50%

2.0m 3.0m

40%

4.0m 5.0m

30%

6.0m 6.9m

20%

10%

1

0% 0

0 0

0.05

0.1

0.15

0.2

0.25

0.3

0.35

0.05

0.1

0.15

0.2

Radial distance [m]

Time-averaged axial H2 flux [kg/m2.s] Fig. 11. The time-averaged axial profile of H2 flux in the optimum case (run 33).

Fig. 12. The time-averaged radial profiles of H2 purity at different heights in the optimum case (run 33).

1055

K. Phuakpunk et al. / Chemical Engineering Science 192 (2018) 1041–1057

and quite uniform distribution and accompany with the axial profile of volume fraction of the total solid discussed previously in Fig. 13, thus the fluidization occurred in the system of run 33 was in fast fluidization regime. In Fig. 16, the time-averaged radial profile of volumetric catalyst to sorbent ratio at different heights were plotted comparing with the ratio where the solid entranced (1.778 vol./vol.). The results indicated that profiles near the left wall had little higher ratio than the inlet ratio but the ratios near the right wall were little less than the inlet ratio. This is because solids entranced from the channel on the right wall with normal direction to the wall so the heavier catalyst would flow directly to the left wall more than the sorbent. While the lighter sorbent would be lifted up at right zone more than the catalyst. However, the ratios had only slightly change

7

6 Catalyst Sorbent

5

Height [m]

Solid

4

3

2

0 0

0.05

0.1

0.15

0.2

0.25

0.3

0.35

0.4

Time- and area-averaged axial volume fracon [-] Fig. 13. The time- and area-averaged axial profiles of volume fraction of solid phases in the optimum case (run 33).

Time-averaged radial volume fracon [-]

0.1 m height 1

1.0 m height

0.6

0.6

0.5

0.5

0.4

0.4

0.3

0.3

0.2

0.2

0.1

0.1

Catalyst Sorbent Solid

0

0 0

0.05

0.1

0.15

0

0.2

0.1

0.15

0.2

3.0 m height

2.0 m height Time-averaged radial volume fracon [-]

Fig. 15 shows the time-averaged radial profile of volume fraction of each solid phase at different heights. At the height over the solid entrance (0.1 m height), the total solid were dense near the solid entrance on the right but dilute at distance near the left wall. When the solids were flowed up to 4.0 m height, the volume fraction of solids at the right wall would gradually decrease to balance with at the left wall. In this height range (0.1–4.0 m height), the solids near the both wall were obviously denser than in the middle. These profiles were core-shell formation of fast fluidization. At heights over 5.0 m, the radial profiles had very low fraction

0.05

Radial distance [m]

Radial distance [m]

0.6

0.6

0.5

0.5

0.4

0.4

0.3

0.3

0.2

0.2

0.1

0.1

Catalyst Sorbent Solid

0

0 0

0.05

0.1

0.15

0

0.2

0.05

0.1

0.15

0.2

Radial distance [m]

Radial distance [m] 7

Time-averaged radial volume fracon [-]

4.0 m height 6 axial inlet

Height [m]

5

4

5.0 m height

0.6

0.6

0.5

0.5

0.4

0.4

0.3

0.3

0.2

0.2

0.1

0.1

Catalyst Sorbent

0

Solid

0 0

0.05

0.1

0.15

0.2

0

Radial distance [m]

3

0 1.6

1.65

1.7

1.75

1.8

1.85

1.9

1.95

2

Time- and area-averaged axial catalyst to sorbent rao [vol./vol.] Fig. 14. The time- and area-averaged axial profile of volumetric catalyst to sorbent ratio in the optimum case (run 33).

Time-averaged radial volume fracon [-]

1

0.1

0.15

0.2

6.9 m height

6.0 m height 2

0.05

Radial distance [m]

0.6

0.6

0.5

0.5

0.4

0.4

0.3

0.3

0.2

0.2

0.1

0.1

Catalyst Sorbent Solid

0

0 0

0.05

0.1

0.15

Radial distance [m]

0.2

0

0.05

0.1

0.15

0.2

Radial distance [m]

Fig. 15. The time-averaged radial profiles of volume fraction of solid phases at different heights in the optimum case (run 33).

1056

K. Phuakpunk et al. / Chemical Engineering Science 192 (2018) 1041–1057

The hydrodynamics of the optimum case indicated that SESMR could be nearly completed since 5.0 m height thus the 7.0 m height of the riser was sufficient. That was because the axial distribution of solid was sufficiently dense below the 5.0 m height without serious segregation of the catalyst and the sorbent while the radial distributions of solid were developed well and quite uniform without serious segregation of the catalyst and the sorbent. SESMR could perform well at almost all zones inside the riser.

2.4 0.1m

Time-averaged radial catalyst to sorbent rao [vol./vol.]

2.3

1.0m 2.0m

2.2

3.0m 4.0m

2.1

5.0m 6.0m

2

6.9m

Acknowledgements

inlet

1.9

This study was funded by the Ratchadapisek Sompoch Endowment Fund (2016), Chulalongkorn University (CU-59-003-IC). During the research, this study had been also supported by the Thailand Research Fund (RSA5980052); the Center of Excellence in Catalysis and Catalytic Reaction Engineering, Chulalongkorn University; the Center of Excellence on Petrochemical and Materials Technology, Chulalongkorn University; and the Department of Mechanical Engineering, Chulalongkorn University

1.8 1.7 1.6 1.5

References

1.4

0

0.05

0.1

0.15

0.2

Radial distance [m] Fig. 16. The time-averaged radial profiles of volumetric catalyst to sorbent ratio at different heights in the optimum case (run 33).

(were in 1.6–2.1 vol./vol.). Therefore, the segregation between the catalyst and the sorbent in radial direction did not seriously affect SESMR and accompany with the axial profile of volumetric catalyst to sorbent ratio discussed previously in Fig. 14, thus SESMR could perform well at almost zone inside the riser and this maximized the H2 flux and made the H2 purity reach the equilibrium.

4. Conclusion In this study, CFD simulation had been used to determine a pilot-scale riser suitable for the highest hydrogen production via SESMR in CFBR system. The Euler-Euler models with the KTGF theory were applied with heterogeneous kinetics of the catalyst and the sorbent in 2D and transient simulation. The 25 full factorial design and ANOVA method were used to investigate effects of design parameters, including the riser diameter, the gas velocity and the solid flux, and reaction parameters, including the inlet temperature and the catalyst to sorbent ratio, on the H2 flux and the H2 purity. The results of 32 runs showed that the SESMR was feasible to operate in the CFB riser because the H2 purity in many cases could reach 99% as theoretical equilibrium and there was no CaO conversion in any case that was over 28% of circulating limitation. There were two similarities of ANOVA results of both the H2 flux and the H2 purity. Firstly, the design parameters were obviously more effective than the reaction parameters. Lastly, the gas velocity, the riser diameter and the solid flux were the most three significant in descending order. Most of main effects/interactions positively affected both the H2 flux and the H2 purity, except the gas velocity and its interactions had negative effects to the H2 purity. To investigate the most suitable system in the studied ranges of the parameters, the possibly highest H2 purity was simulated at 98.58% in dry basis with maximum of 0.301 kg/m2 s H2 flux in the system consisted of the riser with 0.2 m diameter, the inlet temperature of 581 °C, the catalyst to sorbent ratio of 2.54 kg/kg, 200 kg/m2 s of the solid flux and the gas velocity at 6 m/s.

Abanades, J.C., Anthony, E.J., Lu, D.Y., Salvador, C., Alvarez, D., 2004. Capture of CO2 from combustion gases in a fluidized bed of CaO. AICHE J. 50, 1614–1622. Aceves Olivas, D.Y., Baray Guerrero, M.R., Escobedo Bretado, M.A., Marques da Silva Paula, M., Salinas Gutiérrez, J., Guzmán Velderrain, V., López Ortiz, A., CollinsMartínez, V., 2014. Enhanced ethanol steam reforming by CO2 absorption using CaO, CaO*MgO or Na2ZrO3. Int. J. Hydrog. Energy 39, 16595–16607. ANSYS Inc., ANSYS Fluent Theory Guide 15.0, SAS IP Inc., USA, 2013. Arstad, B., Prostak, J., Blom, R., 2012. Continuous hydrogen production by sorption enhanced steam methane reforming (SE-SMR) in a circulating fluidized bed reactor: Sorbent to catalyst ratio dependencies. Chem. Eng. J. 189–190, 413– 421. Barelli, L., Bidini, G., Gallorini, F., Servili, S., 2008. Hydrogen production through sorption-enhanced steam methane reforming and membrane technology: A review. Energy 33, 554–570. Chalermsinsuwan, B., Samruamphianskun, T., Piumsomboon, P., 2014. Effect of operating parameters inside circulating fluidized bed reactor riser with ring baffles using CFD simulation and experimental design analysis. Chem. Eng. Res. Des. 92, 2479–2492. Chao, Z., Wang, Y., Jakobsen, J.P., Fernandino, M., Jakobsen, H.A., 2012. Numerical investigation of the sorption enhanced steam methane reforming in a fluidized bed reactor. Energy Procedia 26, 15–21. Comas, J., Laborde, M., Amadeo, N., 2004. Thermodynamic analysis of hydrogen production from ethanol using CaO as a CO2 sorbent. J. Power Sources 138, 61– 67. Cotton, A., Patchigolla, K., Oakey, J.E., 2013. Overview of, and experimental methodology for sorption enhanced hydrogen production. Energy Procedia 37, 2232–2244. Di Carlo, A., Bocci, E., Zuccari, F., Dell’Era, A., 2010. Numerical investigation of sorption enhanced steam methane reforming process using computational fluid dynamics Eulerian-Eulerian code. Ind. Eng. Chem. Res. 49, 1561–1576. Gidaspow, D., 1994. Multiphase Flow and Fluidization: Continuum and Kinetic Theory Description. Academic Press, UK. Harrison, D.P., 2008. Sorption-enhanced hydrogen production: A review. Ind. Eng. Chem. Res. 47, 6486–6501. Herce, C., Cortés, C., Stendardo, S., 2017. Computationally efficient CFD model for scale-up of bubbling fluidized bed reactors applied to sorption-enhanced steam methane reforming. Fuel Process. Technol. 167, 747–761. Joensen, F., Rostrup-Nielsen, J.R., 2002. Conversion of hydrocarbons and alcohols for fuel cells. J. Power Sources 105, 195–201. Johnsen, K., Grace, J.R., Elnashaie, S.S.E.H., Kolbeinsen, L., Eriksen, D., 2006. Modeling of sorption-enhanced steam reforming in a dual fluidized bubbling bed reactor. Ind. Eng. Chem. Res. 45, 4133–4144. Johnsen, K., Ryu, H.J., Grace, J.R., Lim, C.J., 2006. Sorption-enhanced steam reforming of methane in a fluidized bed reactor with dolomite as CO2-acceptor. Chem. Eng. Sci. 61, 1195–1202. Koumpouras, G.C., Alpay, E., Lapkin, A., Ding, Y., Stepanek, F., 2007. The effect of adsorbent characteristics on the performance of a continuous sorptionenhanced steam methane reforming process. Chem. Eng. Sci. 62, 5632–5637. Kunii, D., Levenspiel, O., 1991. Fluidization Engineering. Butterworth-Heinemann, USA. Kunii, D., Levenspiel, O., 1997. Circulating fluidized-bed reactors. Chem. Eng. Sci. 52, 2471–2482. Lin, J.S., Chen, M.M., Chao, B.T., 1985. A novel radioactive particle tracking facility for measurement of solids motion in gas fluidized beds. AICHE J. 31, 465–473. Montgomery, D.C., 2012. Design and Analysis of Experiments. John Wiley and Sons, USA.

K. Phuakpunk et al. / Chemical Engineering Science 192 (2018) 1041–1057 Ochoa-Fernandez, E., Rusten, H.K., Jakobsen, H.A., Ronning, M., Holmen, A., Chen, D., 2005. Sorption enhanced hydrogen production by steam methane reforming using Li2ZrO3 as sorbent: Sorption kinetics and reactor simulation. Catal. Today 106, 41–46. Prajongkan, Y., 2011. Three-dimensional simulation of hydrodynamics in riser of circulating fluidized bed reactor, Thesis of the degree of Master of Science, Chulalongkorn University. Ranade, V.V., 2002. Computational Flow Modeling for Chemical Reactor Engineering. Academic Press, UK. Rodríguez, N., Alonso, M., Abanades, J.C., 2011. Experimental investigation of a circulating fluidized-bed reactor to capture CO2 with CaO. AICHE J. 57, 1356– 1366. Samruamphianskun, T., Piumsomboon, P., Chalermsinsuwan, B., 2012. Effect of ring baffle configurations in a circulating fluidized bed riser using CFD simulation and experimental design analysis. Chem. Eng. J. 210, 237–251. Sánchez, R.A., Jakobsen, H.A., 2012. Simulation of sorption enhanced steam methane reforming and chemical looping reforming in a circulating fluidized bed reactor. In: Ninth International Conference on CFD in the Minerals and Process Industries. Melbourne, Australia, CSIRO, pp. 1–6. Sánchez, R.A., Solsvik, J., Jakobsen, H.A., 2012. Modeling and simulation of cold flow fluidized bed reactors. Energy Procedia 26, 22–30. Sánchez, R.A., Chao, Z., Solsvik, J., Jakobsen, H.A., 2012. One dimensional two-fluid model simulations of the SE-SMR process operated in a circulating fluidized bed reactor. Procedia Eng. 42, 1282–1291. Sánchez, R.A., Chao, Z., Solsvik, J., Jakobsen, H.A., 2013. An investigation of the heat integration between the two riser units constituting a circulating fluidized bed reactor for the SE-SMR process. Energy Procedia 37, 1218–1227. Solsvik, J., Chao, Z., Jakobsen, H.A., 2012. A one-dimensional two-fluid gas–solid model applied to fluidized bed reactors: The SMR and SE-SMR processes. Procedia Eng. 42, 283–294.

1057

Solsvik, J., Chao, Z., Sánchez, R.A., Jakobsen, H.A., 2014. Numerical investigation of steam methane reforming with CO2-capture in bubbling fluidized bed reactors. Fuel Process. Technol. 125, 290–300. Sun, P., Grace, J.R., Lim, C.J., Anthony, E.J., 2008. Determination of intrinsic rate constants of the CaO–CO2 reaction. Chem. Eng. Sci. 63, 47–56. Wang, Y., Chao, Z., Jakobsen, H.A., 2010. CFD modelling of CO2 capture in the SESMR process in the fluidized bed reactors. Chem. Eng. Trans. 21, 601–606. Wang, Y., Chao, Z., Jakobsen, H.A., 2010. 3D Simulation of bubbling fluidized bed reactors for sorption enhanced steam methane reforming processes. J. Nat. Gas Sci. Eng. 2, 105–113. Wang, Y., Chao, Z., Chen, D., Jakobsen, H.A., 2011. SE-SMR process performance in CFB reactors: Simulation of the CO2 adsorption/desorption processes with CaO based sorbents. Int. J. Greenh. Gas Control 5, 489–497. Wang, J., Wang, Y., Jakobsen, H.A., 2014. The modeling of circulating fluidized bed reactors for SE-SMR process and sorbent regeneration. Chem. Eng. Sci. 108, 57– 65. Wu, Y.J., Li, P., Yu, J.G., Cunha, A.F., Rodrigues, A.E., 2013. Sorption-enhanced steam reforming of ethanol on NiMgAl multifunctional materials: Experimental and numerical investigation. Chem. Eng. J. 231, 36–48. Wu, Y.J., Li, P., Yu, J.G., Cunha, A.F., Rodrigues, A.E., 2014. Sorption-enhanced steam reforming of ethanol for continuous high-purity hydrogen production: 2D adsorptive reactor dynamics and process design. Chem. Eng. Sci. 118, 83–93. Xiu, G., Li, P., Rodrigues, A.E., 2002. Sorption-enhanced reaction process with reactive regeneration. Chem. Eng. Sci. 57, 3893–3908. Xu, J., Froment, G.F., 1989. Methane steam reforming, methanation and water-gas shift: I. Intrinsic kinetics. AICHE J. 35, 88–96. Yeoh, G.H., Tu, J., 2010. Computational Techniques for Multi-phase Flows. Elsevier Ltd, UK.