Generalized model development for a cryo-adsorber and 1-D results for the isobaric refueling period

international journal of hydrogen energy 35 (2010) 3598–3609

Available at www.sciencedirect.com

journal homepage: www.elsevier.com/locate/he

Generalized model development for a cryo-adsorber and 1-D results for the isobaric refueling period V. Senthil Kumar a,*, Sudarshan Kumar b a

India Science Lab, General Motors Global R & D, Creator Building, International Technology Park, Bangalore 560066, India Chemical Sciences and Material Systems Lab, General Motors Global R & D, Warren Technical Center Campus, 30500 Mound Road, Warren, MI 48090, USA b

article info

abstract

Article history:

We have developed 3-D model equations for a cryo-adsorption hydrogen storage tank,

Received 1 December 2009

where the energy balance accommodates the temperature and pressure variation of all the

Received in revised form

thermodynamic properties. We then reduce the 3-D model to the 1-D isobaric system and

7 January 2010

study the isobaric refueling period, for simplified geometry and charging conditions. The

Accepted 8 January 2010

hydrogen capacity evolution predicted by the 1-D axial bed model is significantly different

Available online 13 February 2010

than that predicted by the lumped-parameter model because of the presence of sharp temperature gradients during refueling. The 1-D model predicts a higher hydrogen capacity

Keywords:

than the lumped-parameter model. This observation can be rationalized by the fact that

Hydrogen storage modeling

a bed with temperature gradients on equilibration should desorb gas, whenever the

Cryogenic adsorption

adsorbed phase entropy is lower than the gas phase entropy. The 1-D analysis of the

Metal-organic frameworks

isobaric refueling period does not show any significant difference in hydrogen capacity

Refueling

evolution among the axial, single and multicartridge annular bed designs. Hence, a multicartridge annular design, though giving a slightly lower pressure drop, does not offer any heat and mass transfer enhancement over the single cartridge design. And, the single cartridge annular design appears to be optimal. ª 2010 Professor T. Nejat Veziroglu. Published by Elsevier Ltd. All rights reserved.

1.

Introduction

On-board storage of hydrogen by adsorption at low temperatures (around 77 K, liquid nitrogen temperature) and moderately high pressures (less than 60 bar) is considered viable and competitive with other storage technologies including liquid hydrogen, compressed gas, and metallic or complex hydrides [1]. Nanoporous metal-organic frameworks (MOFs) as adsorbents offer good gravimetric capacity and fast and reversible kinetics. For example, MOF-5 has a reversible capacity of about 6 wt% at 77 K and 50 bar [2]. MOF-5 as an adsorbent material has been studied extensively and we consider it to assess the tank performance in our previous [3] and current works.

Consider a fuel-tank with an initial operating condition of 20 bar and about 80 K, and a fuel cell operating at about 4 bar. The four processes occurring in a cryo-adsorber fuel-tank are:  Refueling – A depleted fuel-tank at low pressure (say 4 bar or 1.1 bar) and higher temperature (say 120–140 K) is charged with hydrogen. Heat released during adsorption is removed by the recirculation of cool hydrogen gas entering at 80 K or below.  Discharge – When the vehicle is in motion, the fuel cell stack demands are met by desorbing hydrogen from the cryoadsorber bed. Desorption is enhanced by the recirculation of hot hydrogen gas. The tank pressure may eventually drop to say 4 bar or 1.1 bar.

* Corresponding author. Tel.: þ91 80 41984567; fax: þ91 80 41158262. E-mail addresses: [email protected], [email protected] (V.S. Kumar). 0360-3199/$ – see front matter ª 2010 Professor T. Nejat Veziroglu. Published by Elsevier Ltd. All rights reserved. doi:10.1016/j.ijhydene.2010.01.025

international journal of hydrogen energy 35 (2010) 3598–3609

Notation T, P temperature and pressure, K, bar L, R, dp length and radius of the adsorbent bed, pellet diameter, m mass and volume of the adsorbent bed, kg, m3 ms, Vb U, Q superficial gas velocity and volumetric flow rate through the bed, m/s, m3/s 3p, 3b, 3t pellet, bed and total porosities, 3t ¼ 3b þ (1  3b)3p, m3/m3 rs, rp, rb skeletal, pellet and bed densities, rp ¼ rs(1  3p), rb ¼ rp(1  3b), kg/m3 rg, vg gas density and specific volume, vg ¼ 1/rg, kg/m3, m3/kg

 Dormancy & Venting – While the vehicle is parked, due to continuous heat leak into the fuel-tank, its pressure builds up to vent pressure, say 25 bar, and gas is vented to the fuel cell (or atmosphere). The process of venting is also called boil-off, in analogy with liquid hydrogen storage systems. The time period till the start of venting is called dormancy. Interestingly these fuel-tank processes occur over different time scales: refueling over a few minutes, discharge over a few hours, dormancy over a few days, and venting over a few weeks. In [3] we had shown that even the fastest process, i.e. refueling, is quasi-static i.e. local equilibrium conditions prevail. Hence, the slower processes are obviously quasistatic. When the molecular processes are fast, slow processes are expected to have negligible internal gradients and are generally amenable for a lumped-parameter analysis. Hence, we had developed a quasi-static lumped-parameter model for the cryo-adsorber fuel-tank in [3]. However, one might anticipate that higher dimensional models might be required, especially for the fastest process i.e. refueling. The refueling period could be divided into two periods:

3599

mg gas viscosity, Pa s aPg, kTg isobaric thermal expansion coefficient and isothermal compressibility, 1/K, 1/bar Hg, Hq, Hs specific enthalpy of gas, adsorbate and adsorbent, J/kg CPg, CPs specific heat capacity of gas and adsorbent, J/kg/K q, q* excess adsorbate concentration and its equilibrium value, kg H2/kg adsorbent the two Langmuir parameter, kg/kg, 1/bar qm, b heat of adsorption, J/kg H2 adsorbed DHa effective thermal conductivity of the adsorbent leff bed, W/m/K

non-ideal, and its specific heat capacity varies significantly with temperature. Hence, in this work we develop a generic form of the 3-D model, which allows variation in all the thermodynamic properties. We then reduce the 3-D model to the 1-D isobaric system and study the isobaric refueling period, for simplified geometry and charging conditions. We study only the isobaric period of refueling herein, since the isobaric period is much longer than the initial adsorption period. We have solved these 1-D equations using finite element method, using the commercial package COMSOLª.

1. Initial adsorption period. During which a depleted bed at low pressure is exposed to the high pressure supply, there are initial pressure transients due to the hydrodynamic filling of the bed and a rapid initial adsorption leading to a temperature spike in the bed. Typically this period lasts a few seconds. 2. Isobaric adsorption period. During which pressure transients are negligible and the bed is getting cooled by and adsorbs further hydrogen from the recirculating cool hydrogen gas. This period is considered isobaric, since the bed is already filled hydrodynamically in the former period and the pressure drop across the highly porous beds are generally negligible. Typically this period extends over a few minutes.

The model equations for adsorption storage or pressure swing adsorption are widely reported in literature [4,5]. These models generally involve approximations like ideal gas behavior, constant specific heat capacities, constant heat of adsorption etc. However, for example, we know that hydrogen at high pressures and cryogenic temperatures is quite

Fig. 1 – Longitudinal view of an axial packed bed.

3600

international journal of hydrogen energy 35 (2010) 3598–3609

Section 2 describes the three packed-bed designs considered in this work: axial, single cartridge annular and the multicartridge annular beds. Section 3 describes the computation of steady-state pressure drop across these bed designs. Section 4 describes the 3-D model development, and in Section 5 it is reduced to the 1-D isobaric system. Section 6 describes the 1-D axial bed results, and compares these results with the lumped-parameter model predictions. Section 7 describes the results for the annular bed designs and compares the three designs. Section 8 studies the effect of feed flow rate on refueling time. The Appendix describes the various input data used in the 1-D isobaric refueling period simulations. And, the support material details the derivation of the basic forms of transient mass and energy balances used in Section 5.

2.

Packed-bed designs

A fixed mass of adsorbent can be packed within a pressure vessel in various ways. The classical style is the axial bed, shown in Fig. 1, where the gas flows axially through a cylindrical bed. An alternative design, giving lower pressure drop, is the annular backed bed, shown in Fig. 2, where the gas is ideally expected to flow radially. Hence, an annular bed design is often used in industrial adsorption or catalytic systems; see for example the metal hydride beds [6,7], adsorption refrigeration systems [8], catalytic reactors [9], chromatographic columns [10] etc. However, it is well-known, both through simulations and experiments that annular packed beds are more prone to flow maldistribution than axial packed beds [11,12]. The presence and effect of flow maldistribution can be systematically studied only in a two-dimensional model. The CFD studies [13,14] show that having the inlet and outlet ports on the same side, such that the inner and outer flows are antiparallel gives low flow maldistribution. Note that this design is different from that given in Fig. 2. Such design variations will be considered in our future two-dimensional studies. Instead of arranging the whole bed as a single annular bed, one could arrange the adsorbent bed into multiple cartridges, as in Fig. 3. In this work, we study the relative advantages of these three different designs, especially in the context of isobaric refueling period. In particular, we seek to answer the following questions: for a fixed mass of adsorbent, is there any process advantage (in heat and mass transfer sense) in increasing the number of cartridges or is there an optimum number of cartridges? These questions cannot be addressed within the premise of a lumped-parameter model, since all the designs are equivalent at that level, and we need at least a one dimensional model. In our earlier work [3] we had considered a MOF of bulk density 385.05 kg/m3 and a bed of 90.37 kg mass and 0.2347 m3 volume. We consider packed beds of the same mass and volume in this study. For all the designs we fix the bed length at 1 m and back calculate the bed diameter from the bed volume. For single and 7-cartridge annular beds, assuming the distribution diameter to be 1 cm, the outer bed diameters are respectively 0.5467 m and 0.2068 m. While comparing different bed designs, we consider only the hydrogen capacity in the bed, both adsorbed hydrogen and the gaseous hydrogen

Feed

Distribution tube

Packed bed

Collection tube

Insulation

Exit

Fig. 2 – Longitudinal and cross-sectional views of a single annular packed bed.

in its pores. The gaseous hydrogen stored in the dead-spaces around the packed bed but within the pressure vessel (for example, space between the cartridges, header and collector spaces etc.) are not considered.

3. Pressure drops for three packed-bed designs We compute the steady-state pressure drop for a nonadsorbing hydrogen gas flow of 120 g/s through these three beds at 20 bar and 80 K. To compute the pressure drop across

Fig. 3 – Cross-sectional view of a 7-cartridge annular packed bed.

international journal of hydrogen energy 35 (2010) 3598–3609

3601

a bed for an unidirectional axial flow the widely used form of Ergun equation is [15, p. 191],

values for bed density and total porosity, we have the general form of mass balance as:

DP ¼ amg U þ brg U2 ; L

rb

(1)

150ð1  3b Þ2 33b d2p

1:75ð1  3b Þ : 33b dp

(3)

Now, we introduce the quasi-static approximations i.e. the gas and the adsorbate are in local equilibrium. This statement implies that any change in the property of gas or the adsorbate at any location is solely due to the temperature and pressure change at that location. For example, with gas density we have

/

(4)

For a radial flow in a differential element this equation becomes vP  ¼ amg Ur þ brg U2r ; vr

(5)

where Ur ¼

Q 2prL

(6)

is the radial superficial velocity. Note that in a radial flow the superficial velocity varies with radial location. Integrating this equation for an annular bed from the inner radius ri to the outer radius ro gives  Po ¼ Pi  amg

     2  Q ro Q 1 1  brg :  ln ri 2pL 2pL ri ro

(7)

Using this equation, the pressure drop across a single and seven cartridge annular beds are 337.72 Pa and 12.05 Pa, respectively. These calculations show that there is significant reduction in pressure drop as we move from axial bed design to single cartridge annular bed design to seven cartridge annular bed design.

r ; tÞ ¼ rg ½Tð! r ; tÞ; Pð! r ; tÞ; rg ð!

vrg vT vP ¼ rg aPg þ rg kTg ; vt vt vt

/

/

(13)

Similarly the adsorbate concentration at any location is approximated to the equilibrium adsorbate concentration at that location i.e. r ; tÞ; Pð! r ; tÞ; qð! r ; tÞzq ½Tð!

(14)

and vq vq z : vt vt

(15)

The Langmuir adsorption isotherm is q bP ¼ : qm 1 þ bP

(16)

Temperature dependence of the two Langmuir parameters are represented as in [3] as   B b ¼ b0 exp T

(17)

and qm0 ; f ðTÞ

(18)

(19)

with

4.1.

Transient mass balance

f ðTÞ ¼ 1 þ AT2 :

The adsorbate diffusion in and out of the element is generally negligible compared to the convective component, and hence it is typically neglected in the literature [4]. Assuming constant

/

V rg ¼ rg aPg V T þ rg kTg V P:

Development of 3-D model equations

(8)

(12)

and

4.

The basic form of the transient mass balance (derivation details are in the supporting material) is

(11)

and then the temporal and spatial derivatives of density could be written in terms of those of temperature and pressure using the chain rule as:

qm ¼

 / !  v r q þ 3t rg þ V $ U rg ¼ 0: vt b

(10)

rb

As in [3] we assume a bed porosity of 3b ¼ 0.36 and a pellet diameter of dp ¼ 1.0 mm. From the NIST web book [16], we get hydrogen gas density and viscosity at 80 K and 20 bar as rg ¼ 6.2070 kg/m3 and mg ¼ 3.7453  106 Pa s. Using these values, we get a pressure drop of 1417.33 Pa across the axial bed. Ergun equation is generalized for a three dimensional flow in a differential element as [17, p. 48], ! !!  V P ¼ amg U þ brg j U j U :

/ ! vrg ! / vq þ U $ V rg þ rg ðV $ U Þ ¼ 0: þ 3t vt vt

(2)

and b¼

(9)

Expanding the gradient term we have

where a¼

vrg / !  vq þ V $ U rg ¼ 0: þ 3t vt vt

Then, we have    vq vT vq vP þ vT P vt vP T vt  0    f ðTÞ B vT q vP ¼ q þ þ ; 2 ð1 þ bPÞP vt f ðTÞ ð1 þ bPÞT vt

vq ¼ vt



(20)

see [3] for the intermediate steps. Using these expressions in the mass balance, we have

3602

a11

international journal of hydrogen energy 35 (2010) 3598–3609

/ ! vT vP ! / ! / þ a12  rg aPg U $ V T þ rg kTg U $ V P þ rg ðV $ U Þ ¼ 0; vt vt

(21)

vHg vHg vHs vq vP ! / þ rb q þ rb þ rb DHa  3t þ rg U $ V Hg vt vt vt vt vt   /  / vrg / !  vq þ V $ U rg þ V $ leff V T þ Hg rb þ 3t vt vt

3t rg

where a11 ¼ rb q

 0  f ðTÞ B  rg aPg 3t þ f ðTÞ ð1 þ bPÞT2

¼ 0:

Noting that the coefficient of gas enthalpy in the above equation is identically zero from the mass balance equation (9), we have

and  rb q þ rg kTg 3t : ð1 þ bPÞP

 a12 ¼

4.2.

(23)

Transient energy balance

The basic form of the transient energy balance (derivation details are in the supporting material, see also [15, p. 335]) is   v vP / ! 3t rg Hg þ rb qHq þ rb Hs  3t þ V $ U rg Hg vt vt /  / þ V $ leff V T ¼ 0:

vH  vHs vq vP ! / g þ rb þ rb DHa  3t þ rg U $ V Hg 3t rg þ rb q vt vt vt vt /  / þ V $ leff V T ¼ 0:

(24)

vrg vHg vHq vq vHs vP ! / þ 3t Hg þ rb q þ rb Hq þ rb  3t þ rg U $ V Hg 3t rg vt vt vt vt vt vt / !  /  / þ Hg V $ U rg þ V $ leff V T ¼ 0: (25) In this form we see that the gas and adsorbate enthalpies appear as coefficients of derivatives, which at the outset could absurdly mean that the numerical evolution of this equation depends on the reference state assumed for the computation of these enthalpies. We will now demonstrate that it is not so. The adsorbate enthalpy is obtained from the definition of heat of adsorption as

vP  vHg vT ¼ CPg þ vg 1  aPg T ; vt vt vt / /  V Hg ¼ CPg V T þ vg 1  aPg T V P;

(32)

/

vHs vT vP ¼ CPs þ vs ; vt vt vt

(33)

(34)

see for example Section 2.1 of Ahluwalia and Peng [18]. Using these results in Eq. (30) and simplifying we have,  ! / vT vP ! / þ a22 þ rg CPg U $ V T þ 1  aPg T U $ V P vt vt /  / þ V $ leff V T ¼ 0; a21

(35)

where  0    f ðTÞ B ; þ a21 ¼ 3t rg þ rb q CPg þ rb CPs  q rb DHa f ðTÞ ð1 þ bPÞT2 (36) and      q rb DHa a22 ¼ 3t þ rb vg q 1  aPg T þ ð1  23t Þ þ : ð1 þ bPÞP

(37)

(26) We have used

r b vs ¼ (27)

Using these results, we have  vq vrg vHg vHg vHs vP þ 3t Hg þ rb q þ rb Hg þ DHa  3t þ rb vt vt vt vt vt vt / !  /  / ! / þ rg U $ V Hg þ Hg V $ U rg þ V $ leff V T

3t rg

¼ 0:

(31)

and hence, we have

This study assumes a constant average heat of adsorption, independent of temperature and pressure. Hence, vHq vHg z : vt vt

(30)

In this form, all the enthalpies appear as derivatives. Hence, the numerical evolution of the above equation is obviously independent of the reference states assumed for computing the enthalpies. Now, as in mass balance, let us introduce the quasi-static approximations. Assuming local equilibrium we have r ; tÞ ¼ Hx ½Tð! r ; tÞ; Pð! r ; tÞ; Hx ð!

Since the gas–solid contact area is very large in porous packed beds, local thermal equilibrium could typically be assumed for physisorption scenarios and hence one single energy balance could be used, instead of a separate energy balance for the gas and solid phases, further details are available in [3, Section 4]. Assuming constant bed density and total porosity, the above equation is expanded as follows:

Hq ¼ Hg þ DHa :

(29)

(22)

(28)

Collecting all the terms having gas enthalpy as their coefficient gives

rb ¼ 1  3t rs

(38)

during the above simplifications. Note that this development of transient energy balance does not contain the heat leak into the system. The heat leak term should be considered in the energy balance for the structural steel domain and not the adsorbent bed domain, with a heat flux continuity boundary condition at the interface allowing the heat exchange between the two domains. The energy balance developed above is only for an element in the adsorbent bed. For the system considered herein, during refueling or discharge, the heat effects due to adsorption or

international journal of hydrogen energy 35 (2010) 3598–3609

desorption are significantly larger than the heat leak into to the tank. Hence, the heat leak can be neglected. For example, for an adsorbent with 4 kJ/mol heat of adsorption, adsorbing 5 kg of hydrogen in 5 min, the average rate of heat generation is about 33.3 kW. While typical heat leaks are of the order of a few Watts. Thus, for these two fast processes, i.e. refueling and discharge, it is not necessary to draw two separate energy balances for the structural steel and adsorbent bed domains. However, to estimate the impact of heat leak, one could ‘‘numerically’’ add that heat to the feed gas and study the bed dynamics. For the slow processes, i.e. dormancy and venting, the lumped-parameter description uses a single energy balance which contains the heat leak term. In fact in these processes, the heat leak drives the tank dynamics.

4.3.

Steady-state momentum balance

General solution procedure

4.4.1.

Three-variable formulation

We solve the mass and energy balances, Eqs. (21) and (35), ! coupled with Ergun equation (4) for the variables ð U ; T; PÞ. Simplistically, it could be stated that the mass balance describes the flow rate or velocity changes due to adsorption, the energy balance describes the temperature changes due to the heat released on adsorption, and the Ergun equation describing the pressure drop due to flow. The coupled effects like pressure changes in the gas phase due to adsorption, flow rate or density changes due to temperature gradient etc. are of course built into the model exactly.

4.4.2.

Two-variable formulation

Ergun equation expresses pressure gradient in terms of the superficial velocity. It can be inverted to express superficial velocity in terms of the pressure gradient as follows. Taking modulus of Eq. (4) gives ! ! jV Pj ¼ amg j U j þ brg j U j2 : /

(39)

! This quadratic equation in j U j is solved to get ! j U j ¼ amg þ

/

þ 4brg jV Pj

2brg :

(40)

Here, we use the positive sign for the square root so that the modulus is non-negative. This expression is substituted in the vector equation and rearranged to get / ! !! ! ! V P ¼ amg U þ brg j U j U ¼ amg þ brg j U j U rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi / 1 0 amg þ a2 m2g þ 4brg jV Pj ! @ AU: ¼ 2



Hence, we have

a2 m2g þ 4brg jV Pj:

This expression is used in the mass and energy balances Eqs. (21) and (35), to eliminate superficial velocity in terms of pressure gradient. Then, we solve the coupled mass and energy balances for the variables ( T,P).

4.5.

Comparison with literature

The mass balance and steady-state momentum equations described in Mota et al. [4] are identical to our Eqs. (9) and (4). The energy balance equation of [4], in our notation is  i /  v h vq vP vT ! 3t rg þ rb q T þ rb DHa  3t þ rb CPs þ CPg V $ Trg U vt vt vt vt /  /  V $ leff V T ¼ 0: (43) We will recover the above equation from Eq. (24) by stating the assumptions systematically. Using the adsorbate enthalpy definition Eq. (26), assuming a constant average heat of adsorption Eq. (27), and a constant bed density, Eq. (24) simplifies to Hg

  vH v vq vHs vP g þ rb DHa þ rb  3t 3t rg þ rb q þ 3t rg þ rb q vt vt vt vt vt  /  / ! / (44) þ V $ U rg Hg  V $ leff V T ¼ 0:

Assuming that the effect of pressure on the gas and solid enthalpies are negligible in the given operating range of pressures, and that the specific heat capacities of the gas and solid are constants i.e. they do not vary with temperature, we have, dHg ¼ CPg dT;

(45)

integrating it from a reference temperature and enthalpy gives Hg ¼ Hg0 þ CPg ðT  T0 Þ:

(46)

Similarly, for the adsorbent enthalpy we have dHs ¼ CPs dT:

(47)

Using these results in Eq. (44), and simplifying we have

rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi a2 m2g

(42)

/

amg þ

CPg

Using Ergun equation (4) i.e. the steady-state momentum equation locally for the pressurization of adsorber gives numerical results comparable to that of the unsteady mechanical energy balance (which is merely a dot product of the unsteady momentum balance equation and velocity) [19]. Hence, in adsorption literature it is common to use the Ergun equation locally even when flow and pressure transients exist.

4.4.

/ ! U ¼ 2ð V PÞ rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi

3603

ð41Þ

  / !   i  v v h Hg0 CPg T0 3t rg þrb q þ V $ U rg þCPg 3t rg þrb q T vt vt  / ! vq vT vP þrb DHa þrb CPs 3t þCPg V $ U rg T vt vt vt /  / (48)  V $ leff V T ¼0: Note that the above equation is identical to that in [4] except for the first term, which is identically zero from the mass balance equation (9). Thus, our model reduces to that of [4] when the effect of pressure on enthalpies is negligible and specific heat capacities and heat of adsorption are constant. The model developed herein is quite generic with rigorous accounting of the temperature and pressure dependence of all the thermodynamic properties.

3604

international journal of hydrogen energy 35 (2010) 3598–3609

5. Reduction of 3-D model to 1-D isobaric system For an isobaric system, we drop the Ergun equation (4) and the temporal and spatial derivatives of pressure in the mass and energy balances (21) and (35), and solve the coupled mass and ! energy balances for ð U ; TÞ i.e. mass balance describes the flow rate or velocity changes due to adsorption and energy balance describes the temperature changes due to the heat released on adsorption. For the isobaric system, mass balance equation (21) reduces to a11

/ ! vT ! /  rg aPg U $ V T þ rg ðV $ U Þ ¼ 0; vt

(49)

where a11 is given by Eq. (22). For unidirectional axial flow, it reduces to a11

vT vT vUx ¼ 0:  rg aPg Ux þ rg vx vt vx

(50)

More complicated boundary conditions have been proposed in the literature; see for example [21,22]. However, we observe that the Danckwerts’ boundary condition is used widely in the adsorption literature, see for example [23,24].

6.

Results & discussion for axial bed design

The whole bed is initially at 140 K. As the 80 K feed gas starts blowing through the bed simultaneous cooling and adsorption takes place. For example, in the temperature evolution in Fig. 4 at t ¼ 0.5 min, one could see three regions existing within the bed. The region 0–A is the part of the bed which has already cooled to 80 K and it has got saturated to its capacity at 80 K. The region A–B is the region where simultaneous cooling and adsorption are currently taking place. Here, the temperatures are between the initial bed temperature and the feed gas temperature. As the gas flows through this region it gets heated up. At about point B, the gas gets heated up to 140 K,

For unidirectional radial flow, in cylindrical coordinates, it reduces to a11

vT vT 1 v  rg aPg Ur þ rg ðrUr Þ ¼ 0: vt vr r vr

a (51)

140 t = 0.5 min

B

Expanding the last term gives 1.0 min

120

(52)

The energy balance equation (35), with the additional assumption of constant effective thermal conductivity, reduces to

1.5 min

T (K)

vT vT vUr rg Ur þ ¼ 0: a11  rg aPg Ur þ rg vr r vt vr

2.0 min

100 2.5 min

a21

vT ! / þ rg CPg U $ V T  leff V2 T ¼ 0; vt

(53)

3.0 min

A 80

where a21 is given by Eq. (36). For unidirectional axial flow it reduces to

0

0.5 x (m)

2

a21

vT vT v T þ rg CPg Ux  leff 2 ¼ 0: vt vx vx

(54)

b

120

A

For unidirectional radial flow, in cylindrical coordinates, it reduces to (55)

Expanding the last term and rearranging we have   vT leff vT v2 T a21 þ rg CPg Ur   leff 2 ¼ 0: r vr vt vr

(56)

We solve the coupled mass and energy balance equations in COMSOLª, with the initial conditions: U ¼ 0 and T ¼ T0, all over the domain at t ¼ 0, and the boundary conditions U ¼ Uf and T ¼ Tf at the feed side, and vT=vx ¼ 0 or vT=vr ¼ 0 at the exit side. This Neumann-type exit boundary condition for flow systems is popularly known as the Danckwerts’ boundary condition [20]. Theoretically, this boundary condition can be applied only for a steady system. For an unsteady system, this boundary condition is believed to suppress the gradients near the exit or induce numerical instability. However, in Section 6, we observe that there is no perceptible gradient suppression near exit in our results, even when the adsorption front is moving out of the bed.

3.0 min 2.5 min

mass flow rate (g/s)

  vT vT 1 v vT a21 þ rg CPg Ur  leff r ¼ 0: vt vr r vr vr

1

2.0 min

110

1.5 min 1.0 min

100

t = 0.5 min B 90

0

0.5 x (m)

1

Fig. 4 – (a) Temperature and (b) mass flow rate evolution in an axial flow bed. At t [ 0.5 min, the region 0–A is the saturated region at feed temperature (80 K), the region A–B is the adsorption zone where the bed is simultaneously cooling (from 140 K to 80 K) and adsorbing, and the region B–1 is the non-adsorbing region at initial bed temperature (140 K).

3605

international journal of hydrogen energy 35 (2010) 3598–3609

the initial bed temperature. Hence, beyond point B, the gas just passes through the bed without cooling the bed any further. Thus, in the region B–1 the bed is saturated at its initial temperature. As time progresses region A–B, the adsorption zone, moves downstream, leaving the bed upstream fully saturated at 80 K. In the flow rate evolution within the bed in Fig. 4, the same inferences could be seen. In the adsorption zone A–B flow rate decreases due to adsorption, and in the saturated zones 0–A and B–1 the flow rate remains constant. Once the bed is fully saturated, gas merely passes through the bed without change in flow rate. Fig. 5 shows the gas density and velocity evolution within the bed. Density changes almost by a factor of two within the bed, showing that the temperature dependence of thermodynamic properties should be properly accounted. In the adsorption zone, it is seen that the gas velocity increases steeply due to sharp rise in temperature and consequent fall in density, occurring within the adsorption zone. Consider a circular slice of the axial bed at x, of elemental thickness Dx. Its total hydrogen content is obtained by summing up the adsorbed phase and gas phase contributions as: pR2 Dxrb q ½Tðx; tÞ; Pðx; tÞ þ pR2 Dx3t rg ½Tðx; tÞ; Pðx; tÞ:

(57)

Once we have the temperature and pressure evolution of the bed, integrating this expression over the entire length of the bed gives the total hydrogen content of the bed at any time. We define the refueling time as the time at which the total hydrogen content of the bed is 5 kg (as mentioned in Section 2, we do not consider the gaseous hydrogen storage in the deadspaces). For the isobaric axial bed, with a 120 g/s feed at 80 K and 20 bar, we found the refueling time to be 135.8 s or 2.26 min. From Fig. 4, we can see that even by the end of refueling significant temperature gradients (about 40 K) exist within the bed. We will return to this observation later in this section.

6.1. Comparison of 1-D axial bed model and 0-D model results Fig. 6 shows the hydrogen capacity evolution predicted by the 1-D axial bed model and by the lumped-parameter model for the isobaric refueling period. The respective isobaric refueling times are 5.86 and 2.26 min. (Note that this 0-D result is without the steel thermal mass and heat leak into the system, so that the 0-D and 1-D results can be compared on a common basis.) Such a large discrepancy in the hydrogen capacity evolution and hence the refueling time is due to the significant temperature gradients existing within the bed during the isobaric refueling period, as seen in Fig. 4. We expect lesser disparity between the 0-D and the 1-D models for slower processes like discharge, dormancy and venting, since the temperature gradients are expected to be milder. In order to check this large discrepancy between the 0-D and 1-D models during the isobaric refueling period, we implemented a tanks-in-series model (not reported here). A 0-D model is a single tank model, while the 1-D model is equivalent to a sufficiently large number of tanks-in-series. In fact, serially coupled 0-D models lead to a tanks-in-series model. A tanks-in-series model showed that the two curves shown in Fig. 6 can be smoothly interpolated as the number of tanks increases from 1 to about 25.

2

H load (kg)

5

4

3

0−D

2

1−D Axial bed 1

0

2

4

6

time (mins)

Fig. 5 – (a) Density and (b) superficial gas velocity evolution in an axial bed.

Fig. 6 – Hydrogen capacity evolution predicted by the 0-D model and 1-D axial bed model during the isobaric refueling period.

3606

international journal of hydrogen energy 35 (2010) 3598–3609

from the bed and back into the gas phase, resulting in a slight increase in the gas pressure. Earlier in this section, we noted that even at the end of refueling significant temperature gradients, up to 40 K temperature difference exists within the bed. If the tank is isolated after refueling, with the temperature gradient slowly equilibrating, as we note in this section, due to desorption back into the gas phase there can be a slight overshoot in the tank pressure.

7. Results for annular bed designs and comparison of the three bed designs Figs. 8 and 9 show the simulation results for the single cartridge annular bed. Similar result hold seven cartridge annular bed (figures not shown), except that one seventh of the total flow is apportioned to each cartridge. The phenomenon of simultaneous cooling and adsorption in annular beds is similar to that in the axial bed, except for the geometric effects. From Fig. 8, it can be seen that the adsorption zone moves rapidly near the inner diameter and slowly near the periphery. In Fig. 9, at long times the superficial gas velocity approaches a rectangular hyperbola, since for the non-adsorbing case the radial velocity at any location is given by Eq. (6).

a

140 t = 0.5 min

120

1.0 min

T (K)

Interestingly, from Fig. 6, we see that a 0-D model always predicts lower hydrogen capacity than a 1-D model. The basic distinction between the 0-D and 1-D models is the presence of temperature gradients in the latter. The presence of temperature gradient in a 1-D model facilitates higher hydrogen capacity than the 0-D model. This can be demonstrated by the heuristic model shown in Fig. 7. Consider a kilogram of MOF at 80 K and another kilogram of MOF at 120 K being brought in thermal contact. Assuming that their specific heat capacities do not change significantly, we could say that the equilibrium temperature of the composite block (neglecting for the moment, the energetics of adsorption or desorption) will be about 100 K. From the adsorption data in the Appendix, one could compute the excess adsorbed hydrogen at 20 bar in 1 kg MOF at 80 K as 49.73 g, 1 kg of MOF at 120 K as 19.15 g and 2 kg of MOF at 100 K as 61.69 g. This calculation indicates that 7.19 g of hydrogen is returned to the gas phase when the temperature gradients equilibrate within the composite block. Note that this heuristic model does not capture the ‘‘secondary’’ effects which occur in reality, but illustrates the fact that a bed with gradients has more adsorbed hydrogen than a bed at its equilibrated temperature, and thus it explains qualitatively the discrepancy between the lumped-parameter and 1-D models. We can understand the above fact from a thermodynamic perspective too. If an adsorbent bed with internal temperature gradients is isolated from its environment, the second law of thermodynamics dictates that, on equilibration the system entropy should increase. In most cases of gas physisorption, the entropy of the adsorbed phase is lesser than that of the gas phase. (If the adsorbed molecule dissociates on the surface, it could favor an increase in entropy. But compared to the gas phase, typically an adsorbed phase has much lower degrees of translational freedom. Hence, even with dissociative adsorption, more often the adsorbed phase entropy is lower than that of the gas phase, leading to the commonly observed exothermic adsorption. However, in special cases of dissociative adsorption, if the entropy of the adsorbed phase is higher than the fluid phase, there can be an endothermic adsorption, see [25] for example.) Hence, on equilibration of the temperature gradient some quantity of gas gets desorbed

1.5 min

100

2.0 min 2.5 min

80

3.0 min

0.05

0.1

0.15 r (m)

0.2

0.25

b 1 kg MOF at 80 K & 20 bar

19.15 g H2

+

3.0 min

120

1 kg MOF at 120 K & 20 bar

mass flow rate (g/s)

49.73 g H2

2 kg MOF at 100 K & 20 bar

Returned to gas phase

61.69 g H2

7.19 g H2

Fig. 7 – Heuristic model to demonstrate that a bed with temperature gradients contains more adsorbed hydrogen than the bed at its average temperature.

2.5 min 2.0 min

110 1.5 min 1.0 min

100 t = 0.5 min

90

0.05

0.1

0.15 r (m)

0.2

0.25

Fig. 8 – (a) Temperature and (b) mass flow rate evolution in a single cartridge annular bed.

3607

international journal of hydrogen energy 35 (2010) 3598–3609

a

3.0 min

6

5

2.5 min

H load (kg)

5

2

1.5 min

g

3

(kg/m )

2.0 min

1.0 min

4

4

3 Axial Annular−1

2

Annular−7 t = 0.5 min

1

3

b

0

0.05

0.1

0.15 r (m)

0.2

1

0.25

2

3

time (mins)

Fig. 10 – Hydrogen capacity evolution among the different packed-bed designs.

0.6 0.035 t = 0.5 min 0.03

U (m/s)

0.4

0.025

1.0 min

0.02 0.015

0.2

1.5 min 2.0 min 2.5 min

0.01 0.1

0

0.05

0.15

0.1

0.15 r (m)

0.2

3.0 min 0.25

0.2

0.25

flow maldistribution. The material cost increases with the number of cartridges due to an increase in the mass of structural steel (increased diameter of pressure vessel and headers, additional distribution tubes, bed restrainers etc.). Since flow is routed parallel in a multicartridge design, a flow short circuit to the outer header in a cartridge will lead to significant flow maldistribution. As the number of cartridges increases, the radial thickness of the bed in a cartridge decreases and the chances of short circuiting or flow maldistribution increase. Thus, a single cartridge annular bed offers the advantages of radial flow system, without increasing the material cost, fabrication complexity and flow maldistribution as in a multicartridge annular bed design.

Fig. 9 – (a) Density and (b) superficial gas velocity evolution in a single cartridge annular bed.

8. Fig. 10 compares the hydrogen capacity among the three different packed-bed designs under consideration. Surprisingly, there is negligible difference between the three designs in terms of hydrogen capacity evolution. This observation may be explained as follows: since the internal surface area within the pellets is very high, the adsorbent is already in intimate contact with the gas. Hence, rearranging the bed geometry does not make much difference in the rate at which adsorption happens. The primary motivation of this study was to answer the following questions: for a fixed mass of adsorbent, is there any process advantage in increasing the number of cartridges in an annular bed design and is there an optimum number of cartridges during the isobaric refueling period? Since we do not find significant difference in the hydrogen capacity evolution, we infer that there is no significant heat and mass transfer advantage in moving from axial bed to single cartridge annular bed or to seven cartridge annular bed designs. Hence, the only difference in these design enhancements is the decrease in pressure drop, as we noted in Section 3. Increase in the number of cartridges in an annular design results in higher material cost, fabrication complexity and

Refueling times at various feed flow rates

The results presented in above sections are for a constant feed flow rate of 120 g/s. Now we present the refueling time at different feed flow rates, for the three designs under consideration, in Table 1. Especially, in the axial bed results, we see that doubling the feed flow rate approximately halves the refueling time. This scaling can be reasoned as follows: in this isobaric model, effective thermal conductivity of the bed is the only transport property which contains the dimension of

Table 1 – Refueling time for the three bed designs at various feed flow rates. Feed flow rate (g/s)

240 120 60 30 15

Refueling time (s) Axial bed

Single cartridge annular bed

Seven cartridge annular bed

68.0 136.1 272.4 545.6 1093.6

68.2 136.9 275.7 557.8 1140.1

69.5 141.8 294.0 625.0 1380.5

3608

international journal of hydrogen energy 35 (2010) 3598–3609

time, and all the remaining properties are thermodynamic or material properties which do not contain the dimension of time. When the effective thermal conductivity is low, as in MOF, it does not significantly contribute to the thermal diffusion of the adsorption front through the bed, leaving superficial gas velocity or flow rate as the main variable scaling the refueling time. This scaling, though approximate, helps us in estimating the refueling time at arbitrary feed flow rates given the refueling time at a flow rate, since the product of refueling time and feed flow rate is approximately a constant, for a given feed temperature and pressure. Compared to the axial bed, the single cartridge annular bed design has a lower flow path length. As the number of cartridges increases, the radial thickness of the bed and hence the flow path length decreases. As flow path length decreases, the effect of thermal conductivity of the adsorbent on the thermal diffusion of the adsorption front increases. Hence, we can see that the above scaling is less accurate as we move from the axial bed to single cartridge annular bed to seven cartridge annular bed designs.

For all the bed designs, we find that the product of refueling time and feed flow rate is approximately constant for a given feed temperature and pressure, i.e. doubling the feed flow rate halves the refueling time. This scaling, though approximate, helps us in estimating the refueling time at arbitrary flow rates given the refueling time at a flow rate.

Acknowledgments The authors would like to thank Raghunathan K, Michael Herrmann, Scott Jorgensen, Ulrich Eberle, Dieter Hasenauer, Rainer Immel, Gregory Meisner, Anne Dailly and Eric Poirier for technical data and/or valuable suggestions during the course of model development. Senthil Kumar would like to acknowledge constructive comments from his group members at India Science Lab during his presentations while the model development was in progress.

Appendix A. Simulation input data A.1.

9.

Conclusions

We have developed generalized 3-D model equations for the cryo-adsorber, accommodating the temperature and pressure variation of all the thermodynamic properties. We have shown that in the limiting case of negligible temperature and pressure dependence of the thermodynamic properties and a constant heat of adsorption the derived energy balance equation reduces to that reported in the literature. We then reduce the 3-D model to the 1-D isobaric system and study the isobaric refueling period, for simplified geometry and charging conditions. We have considered three packed-bed designs: an axial bed, single cartridge annular bed and a seven cartridge annular bed design. The hydrogen capacity evolution predicted by the 1-D axial bed model is significantly different from that predicted by the lumped-parameter model, since sharp temperature gradients exist within the bed during refueling. We observe that the lumped-parameter model predicts a lower hydrogen capacity than the 1-D model. This observation can be explained by the fact that due to the second law of thermodynamics, on equilibration an isolated bed with temperature gradients should desorb gas whenever the adsorbed phase entropy is lower than the gas phase entropy. We expect lesser disparity for the slower processes, since the temperature gradients are expected to be milder. During the isobaric refueling period, we did not find any significant difference in hydrogen capacity evolution among the axial bed, single cartridge annular bed and multicartridge annular bed designs. From this observation we conclude the following: annular bed design gives a lower pressure drop than an axial bed design but it does not enhance the heat and mass transfer processes. A multicartridge annular bed design gives a lower pressure than a single cartridge annular bed design but again increase in the number of cartridges does not enhance the heat and mass transfer processes. Hence a single cartridge annular bed design offers the advantages of the annular bed design without increasing the material cost, fabrication complexity and flow maldistribution.

Gas data

We have employed the isobaric thermodynamic properties of hydrogen gas at 20 bar from the NIST web book [16], in the temperature range 60–300 K, and fitted the expressions shown below. vg ¼ vg0 þ vg1 T þ vg2 T2 þ vg3 T3 þ vg4 T4

(A.1)

vg0 ¼ 0:042693; vg1 ¼ 0:0028869; vg2 ¼ 5:5122e  6; vg3 ¼ 1:6913e  8; vg4 ¼ 1:9513e  11; where vg is the specific volume in m3/kg and T in K. We get the isobaric thermal expansion coefficient as aPg ¼

CPg

  vg1 þ 2vg2 T þ 3vg3 T2 þ 4vg4 T3 1 vvg ¼ ; vg vg vT P     exp fT2 T i2 þ f3 þ f4   100 exp fT2  1    2  2 T 100 100 þf6 ; þ f5 þ f7 100 T T

(A.2)

 2 f2 ¼ f1  h T

(A.3)

f1 ¼ 4:788937e þ 03; f2 ¼ 5:229983e þ 02; f3 ¼ 1:496674e þ 04; f4 ¼ 1:204579e þ 03; f5 ¼ 8:489895e þ 01; f6 ¼ 5:356621e þ 03; f7 ¼ 2:625154e þ 03: where CPg is the isobaric specific heat capacity in J kg1 K1 and T in K.

A.2.

Adsorption/adsorbent data

For the Langmuir adsorption model, Eqs. (16)–(19), the parameB ¼ 332.0158 K, ters used are: b0 ¼ 2.8164  108 Pa1, 2 qm0 ¼ 11.7134  10 kg H2/kg MOF and A ¼ 131.3231  106 K2. As in [3], the remaining adsorbent parameters are DHa ¼ 4:0 kJ=mol ¼ 2  106 J=kg, leff ¼ 0:32 W=m=K, and

international journal of hydrogen energy 35 (2010) 3598–3609

log10 CPs ¼

8 X

n  cn log10 T ;

(A.4)

n¼0

where (c8, ., c0) ¼ (0.447169361516, 5.888038301049, 32.914887850468, 101.724106321425, 189.552581009475, 217.283997486133, 148.061506476350, 51.560061018285, 6.714743948231), CPs is in J kg1 K1 and T in K. The bed density and porosity values are from [3]: rb ¼ 385:05 kg=m3 ; 3b ¼ 0:36; 3t ¼ 0:7837.

A.3.

Initial and boundary conditions

Using the lumped-parameter model reported in [3], for a feed at 20 bar pressure and 80 K temperature, and a flow rate of 120 g/s, we get a temperature 134.25 K at the end of the initial adsorption period. In a 1-D model we can neglect the thermal mass of the structural components and the heat leak into the tank through them, since the focus is only on obtaining the data for bed dynamics. So we recalculated the lumpedparameter model result for the initial adsorption period neglecting the steel thermal mass and heat leak into the tank. A tank initially at 1.1 bar and 120 K had a temperature of 140.87 K at the end of initial adsorption period i.e. when the tank pressure reached 20 bar. Since we have not yet performed the 1-D analysis for the initial adsorption period, we arbitrarily assume the initial condition for the isobaric refueling period as P ¼ 20  105 Pa, T0 ¼ 140 K, since it is representative of typical operating conditions. For the axial, single cartridge annular and seven cartridge annular designs the flow velocities at the feed side are 0.08237 m/s, 0.6154 m/s and 0.0879 m/s. We observed the following simulation artifact with the Lagrange elements: upstream of the adsorption front there was a small dip in the temperature even below the feed temperature, which is surely unphysical. However, with the use of Hermite cubic elements, this artifact vanished.

Appendix B. Supplementary data Supplementary data associated with this article can be found, in the online version, at doi:10.1016/j.ijhydene.2010.01.025.

references

[1] Zhou L. Progress and problems in hydrogen storage methods. Renew Sustain Energy Rev 2005;9:395–408. [2] Dailly A, Vajo JJ, Ahn CC. Saturation of hydrogen sorption in Zn benzenedicarboxylate and Zn naphthalenedicarboxylate. J Phys Chem B Lett 2006;110:1099–101. [3] Kumar VS, Raghunathan K, Kumar S. A lumped-parameter model for a cryo-adsorber hydrogen storage system. Int J Hydrogen Energy 2009;34:5466–75.

3609

[4] Mota JPB, Rodrigues AE, Saatdjian E, Tondeur D. Dynamics of natural gas adsorption storage systems employing activated carbon. Carbon 1997;35:1259–70. [5] Arumugam BK, Banks JF, Wankat PC. Pressure effects in adsorption systems. Adsorption 1999;5:261–78. [6] Gopal MR, Murthy SS. Prediction of heat and mass transfer in annular cylindrical metal hydride beds. Int J Hydrogen Energy 1992;17:795–805. [7] Fukada S, Morimitsu S, Shimoozaki N. Heat and mass transfer in concentric-annular-tube bed packed with ZrV1. 9Fe0.1 particles. J Alloys Compd 2004;375:305–12. [8] Liu Y, Leong KC. The effect of operating conditions on performance of zeolite/water adsorption cooling systems. Appl Therm Eng 2005;25:1403–18. [9] Yoo CS, Dixon AG. Modeling and simulation of a mixed-flow reactor for ammonia and methanol synthesis. Chem Eng Sci 1988;43:2859–65. [10] Gu T, Tsai GJ, Tsao GT. A theoretical study of multicomponent radial flow chromatography. Chem Eng Sci 1991;46:1279–88. [11] Ponzi PR, Kaye LA. Effect of flow maldistribution on conversion and selectivity in radial flow fixed-bed reactors. AIChE J 1979;25:100–8. [12] Bolton GT, Hooper CW, Mann R, Stitt EH. Flow distribution and velocity measurement in a radial flow fixed bed reactor using electrical resistance tomography. Chem Eng Sci 2004; 59:1989–97. [13] Heggs PJ, Ellis DI, Ismail MdS. The modeling of fluid-flow distributions in annular packed beds. Gas Separ Purif 1994;8: 257–64. [14] Kareeri AA, Zughbi HD, Al-Ali HH. Simulation of flow distribution in radial flow reactors. Ind Eng Chem Res 2006; 48:2862–74. [15] Bird RB, Stewart WE, Lightfoot EN. Transport phenomena. 2nd ed. NY: John Wiley & Sons; 2002. [16] NIST web book. http://webbook.nist.gov/chemistry/fluid/. [17] Kaviany M. Principles of heat transfer in porous media. 2nd ed. NY: Springer-Verlag; 1995. [18] Ahluwalia RK, Peng JK. Dynamics of cryogenic hydrogen storage in insulated pressure vessels for automotive applications. Int J Hydrogen Energy 2008;33:4622–33. [19] Sereno C, Rodrigues A. Can steady-state momentum equations be used in modelling pressurization of adsorption beds? Gas Separ Purif 1993;7:167–74. [20] Danckwerts PV. Continuous flow systems: distribution of residence times. Chem Eng Sci 1953;2:1. [21] Parulekar SJ, Ramakrishna D. Analysis of axially dispersed systems with general boundary conditions – I. Chem Eng Sci 1984;39:1571–9. [22] Lee TT, Wang FY, Newell RB. On the evaluation of the exit boundary condition for the axial dispersion bioreactor system. Chem Eng Technol 1998;21:901–10. [23] Huang WC, Chou CT. Comparison of radial and axial-flow rapid pressure swing adsorption processes. Ind Eng Chem Res 2003;42:1998–2006. [24] Sankararao B, Gupta SK. Modeling and simulation of fixed bed adsorbers (FBAs) for multi-component gaseous separations. Comput Chem Eng 2007;31:1282–95. [25] Thomas JM. The existence of endothermic adsorption. J Chem Educ 1961;38:138–9.