Accurate simplified dynamic model of a metal hydride tank

.


Hydrogen storage
Reversible fuel cell a b s t r a c t As proton exchange membrane fuel cell technology advances, the need for hydrogen storage intensifies.Metal hydride alloys offer one potential solution.However, for metal hydride tanks to become a viable hydrogen storage option, the dynamic performance of practical tank geometries and configurations must be understood and incorporated into fuel cell system analyses.A dynamic, axially-symmetric, multi-nodal metal hydride tank model has been created in Matlab-Simulink Ò as an initial means of providing insight and analysis capabilities for the dynamic performance of commercially available metal hydride Energy and mass balances are solved in cylindrical polar coordinates for a cylindrically shaped tank.The model tank temperature, heat release, and storage volume have been correlated to an actual metal hydride tank for static and transient absorption and desorption processes.A sensitivity analysis of the model was accomplished to identify governing physics and to identify techniques to lessen the computational burden for ease of use in a larger system model.The sensitivity analysis reveals the basis and justification for model simplifications that are selected.Results show that the detailed and simplified models both well predict observed stand-alone metal hydride tank dynamics, and an example of a reversible fuel cell system model incorporating each tank demonstrates the need for model simplification.

Background and motivation
There is a significant need for advanced hydrogen storage technology to enable the use of proton exchange membrane fuel cells in many applications.Fuel cell automobiles, small portable fuel cell devices, uninterruptible power supply systems [2], and Regenerative Fuel Cells (RFC) will rely on stored hydrogen.The current work is motivated by the desire to develop an RFC system, including metal hydride hydrogen storage that is germane to military auxiliary power unit applications.
Regenerative fuel cells are well-suited to military applications because of the desire to use common fuels.Jet Propulsion 8, or JP-8, is a kerosene-based fuel that is the single battlefield fuel for Department of the Army applications including electric power generators, wheeled and tracked vehicles, aircraft, stoves, and heaters.An immense logistical effort is necessary to move the enormous quantities of JP-8 fuel required for military operations into a forward deployed overseas location.Delivered costs of JP-8 to Army combat platforms have been estimated to be at least $30-$40/gallon for overland transport, and greater than $400/gallon for air delivery [3].
Because of the ''single battlefield fuel'' mandate, fuel cell generators or other clean and efficient hydrogen consuming devices cannot be utilized on the battlefield unless the hydrogen is obtained from reformation of JP-8 or by other means with existing battlefield resources.Consequently, an RFC ''black box'' having electricity as the only input and output, is one potential fuel cell technology available to the military.
The key feature of an RFC system is the ability to independently design the system power capability, energy storage capacity, and recharge rate.Power production is determined by fuel cell size, recharge rate by electrolyzer size, and storage capacity by the amount of fuel storage available.Consequently, a fuel cell/electrolyzer system can theoretically achieve a better energy density than even state-of-the-art chemical batteries [4].This could benefit military applications that currently rely on lead-acid batteries for electrical storage, without violating the single fuel mandate.However, this energy storage potential directly relies on the performance of the hydrogen storage system.
In order to evaluate the performance of a fuel cell rechargeable energy system, to garner insight into the dynamic response characteristics of the system and individual components, and to design superior systems in the future, a dynamic model of a metal hydride tank is developed in a modular Matlab-Simulink framework and verified by experiment.The Simulink framework will allow the hydride model to be easily integrated into a regenerative fuel cell system model.
The physics of this dynamic model follows the work of Mayer et al. [1].The model described in this paper is based on a bulk and geometrically resolved first principles approach that solves the dynamic conservation of mass and conservation of energy equations together with appropriate calculations for heat transfer and chemical reactions, drawing on the substantial metal hydride modeling work in the literature.For example, Sun and Deng [5] used Fortran 77 to study a twodimensional hydride bed.Jemni and Nasrallah [6] studied a two-dimensional, cylindrical LaNi 5 model to demonstrate that convection heat transfer between hydrogen gas and solid metal alloy is negligible.With their numeric model, Gopal and Murthy [7] stressed that hydriding performance is intimately coupled to heat transfer within the metal hydride bed.Jemni et al. [8] published much cited experimental data in 1999 to back their experimental determination of hydride conduction heat transfer coefficients.Aldas et al. [9] used PHOENICS code to show with their three-dimensional model that the additional resolution of a third dimension does not significantly affect model results for determining the mass of hydrogen absorbed.Askri et al. [10] showed, using the conservation of momentum equation, that the small expansion volume normally used in metal hydride tanks is insignificant to the hydriding process.The model developed here is similar to these others in its physics, but accels due to the Simulink environment, in its ability to be easily integrated into complex system models.Integration of hydrogen storage models into larger system models is becoming more important as metal hydride applications become more widespread [11].Additionally, the experimental work reported in the literature generally centers on purpose built, ideal, experimental hydride reaction beds cooled by integrated cooling systems [1,7,8]

Experiment description
The test apparatus used in the experimental portion of this research consisted primarily of a small (approximately 1 kg) Ergenics, Inc. aluminum tank containing Hydralloy C5 metal hydride.To support the testing, an analog pressure gauge, two small manually operated needle valves, a digital temperature indicator with thermocouple, several feet of nylon tubing, a supply of research grade hydrogen gas regulated to 13.6 atm (200 psi), and a water bath were used.This equipment was configured as pictured in Fig. 1.Fig. 2 shows a schematic drawing of the experiment.All of the testing was conducted on a seasoned, previously used hydride tank.This ensured that minimal variation would be observed amongst similar tests.Testing was conducted by first measuring the initial mass of the tank to determine the starting quantity of hydrogen stored using a precision electronic balance.With the hydrogen release valve closed, the inlet valve was then opened to allow hydrogen from a large gas cylinder to enter the tank at 13.6 atm.External tank temperature and internal tank pressure were recorded for the duration of the test.The mass of the tank was again measured at the end of the test to determine the total mass of hydrogen stored.The tank was emptied by closing the inlet valve and opening the hydrogen release valve after each test.
Tests were also conducted with the tank cooled by a circulating water bath, pictured in Fig. 1 and shown schematically in Fig. 3.The procedure was identical to the aircooled test, with the exception that the tank was submerged in a continuously circulating water bath.

Model description
Closely following the work of Mayer et al. [1] a hydride tank model was created to determine the tank temperature, the quantity of energy released, and the mass of hydrogen absorbed within each node of a discretized model when given local operating conditions such as temperature, heat transfer rates, and hydrogen inlet pressure.This was accomplished by simultaneously solving the conservation of mass and conservation of energy equations for the gaseous hydrogen and the solid metal alloy in each of the discrete nodes of the model.

Hydrogen conservation of mass
The mass of hydrogen within a node can change due to advection, or via the sorption reaction with the metal alloy according to: where 3 is the porosity of the metal alloy.The rate of the hydrogen sorption reaction per unit volume, _ m reaction , can be found empirically as described later.Darcy's Law is used to determine the gas velocity, v g based on gas pressure, P, dynamic viscosity of the gas, m, and the metal hydride permeability, K: With Eq. ( 2), Eq. ( 1) can be numerically integrated to find the instantaneous density of gaseous hydrogen in the control volume.

Solid alloy conservation of mass
The volume of solid hydride alloy within a node, V s , is equal to the non-porous volume of the node: The mass of the solid, m s , can only change by absorbing or desorbing hydrogen.The conservation equation is consequently written as:

Hydrogen conservation of energy
The gaseous energy of each node can change by four processes: 1. Conduction with gas in surrounding nodes.2. Convection between gas and the solid alloy.
3. Change in mass due to the sorption process.4. Advection (gas transport) with surrounding nodes.
These four transfer methods are accounted for in the following gaseous energy conservation equation: where h sg is the convection coefficient between hydrogen and the metal alloy, k g is the conduction coefficient between hydrogen in each node, and A is the ratio of surface area to volume between gas and solid within the node.The nodal hydrogen temperature can be determined by numerically integrating over time.

Solid alloy conservation of energy
The change in energy of the solid alloy in each node is due to four processes: 1. Conduction with solid in surrounding nodes.
2. Convection with hydrogen gas within node.
3. Heat release or supply during the sorption reaction.4. Mass change via hydrogen sorption.
The energy conservation equation for the solid within each node is therefore: where the enthalpy of reaction, DH, is expressed in J/kg.The nodal solid temperature can be found by numerically integrating.

Rate of hydrogen reaction
The rate at which hydrogen gas is absorbed or desorbed can be modeled by the following empirical equation: where C a is the reaction rate constant, E a is the activation energy, Ris the ideal gas constant, and r f and r s are the saturated and nodal densities of the solid alloy, respectively.P g is the hydrogen pressure within the node and P vH is the pressure predicted by the van't Hoff equation [12], The values of enthalpy of reaction, DH, and the change in entropy for the reaction, DS, are experimentally determined for the given hydride alloy being used [12].Note that the unit of pressure in Eq. ( 8) is atm.

Heat transfer by convection
As shown in Section 2, the tank was tested under both a free convection condition in air and a forced convection condition in a circulating water bath.The free convection heat transfer coefficient of the vertically oriented hydride cylinder was modeled by representing the sides of the tank as a vertical plate.The average convection coefficient, h, for this cooling scenario can be found from the average Nusselt number, Nu, as: where k is the thermal conductivity of air and L is the length of the tank.For laminar flow on a vertical plane, the average Nusselt number is given by [13]: where Ra is the Rayleigh number and Pr is the Prandtl number.
For the circulating water cooling scenario, the average convection coefficient for forced convection on a cylinder in crossflow can also be found from the average Nusselt number: where D is the diameter of the cylinder.The average Nusselt number is given by [13]: where Re is the Reynolds number.

Baseline model
A baseline model was created using Matlab-Simulink software that implemented all of the above physical relations and was discretized into five cylindrical nodes and five axial divisions as shown in Fig. 4. Hydrogen can only enter or exit the tank at one end in a manner that physically represents the inlet of the commercial hydride tank employed in the experimental work.
Integration blocks in Simulink are the most computationally challenging aspect of the model due to the iterative algorithm needed to converge on a solution.For each node, integration is necessary to resolve the solid temperature, the gas temperature, and the gas mass conservation.This baseline model of 25 total nodes therefore contains 75 coupled integration steps.Increasing the model to 6 Â 6 nodes, or 108 coupled integration steps, led to an impractically slow standalone model that would not be practical in a larger, integrated system model as desired.
Initial execution of the 5 Â 5 node model under air-cooling conditions showed negligible temperature and mass absorption gradients between nodes in the axial direction, and relatively large gradients in the radial direction as shown in Figs. 5 and 6 after 900 s of a filling event.Model execution at other times and other conditions showed similar results.
The axial and radial temperature and mass distributions are a result of the aspect ratio of the particular tank being studied (0.125 m long by 0.0254 m radius, 4.92:1).This can be demonstrated by altering the aspect ratio as shown in Fig. 7 where axial temperature variation is plotted against aspect ratio.
As a result of this observation, the model axial discretization was reduced to just one node, and the number of radial nodes was increased to the practical computational maximum limit of 40 nodes, creating a one-dimensional radial model.This may not be a viable modeling option for all tank shapes and sizes.

Model evaluation by data comparison
Hydride tank model simulations were performed with material properties and physical values matching the actual experiment as closely as possible.The model parameter values, as shown in Table 1, were garnered from material property databases, manufacturer specifications, and published work.Manufacturer and online data presented are for the actual alloy used in testing, Hydralloy C5.Published data is for a similar AB 2 alloy, TiMn 1.5 .
The first simulation mimicked the first physical test in which the tank was filled with pure hydrogen at 13.6 atm

Radial Node Axial Node Grams of H2 Absorbed per m3 of Alloy
Fig. 6 -Tank hydrogen absorption levels for each node in the radial and axial directions.
while the tank was cooled only by natural convection with the surrounding air.As shown in Fig. 8, the model utilizing reference DHand DS values under-predicts the tank temperature.In the experiment, the actual tank absorbed 1.3 g of hydrogen, while the model predicted absorption of just 0.20 g.These are clearly poor results.However, because some of the model parameters were estimated from literature sources that were not directly applicable to this particular metal hydride, disagreement was not unexpected.Fig. 9 shows results for filling the tank while it is cooled by a water bath.Again, the same reference values result in an inaccurate absorption prediction.For this case, the actual tank absorbed 3.1 g of hydrogen while the model predicted 2.1 g, which is not as poor a comparison as the air-cooled case.
As a result of findings presented in Figs. 8 and 9 it was decided that use of the current literature model parameters was not adequate.Instead of using the reference values of Table 1 for the enthalpy and entropy of reaction, these values can be calculated experimentally using the van't Hoff equation for a specific single set of experimental conditions [12].On a plot of ln(P) versus 1/T, DH can be calculated from the slope of the line, and DS is found from the y-intercept.Data collected for the tank used in this testing is presented in Fig. 10 in this manner.
Experimentally determined enthalpy and entropy values are calculated as À28,800 J/mol and À112.2J/mol K, respectively.Using these values, the model accurately reproduces experimental results for the natural convection case, filling with 1.3 g of hydrogen.The temperature versus time plots are shown in Fig. 11.Note that these dynamic model predictions much better approximate the experimental observations.
Model results for the water-cooled tank scenario also show improvement when using the experimentally determined values for enthalpy and entropy of reaction over those from published databases.The model accurately predicted 3.1 g of hydrogen absorption and the dynamics of the tank temperature response are well predicted by the model, as shown in Fig. 12.

Sensitivity analysis to improve model runtime
The 40-node model that employs all of the governing physics well matches experimental data, but, requires approximately 0.15 s of execution time for each real second simulated (using a HP Pavillion laptop computer with an Intel Centrino processor running Matlab 7.0, Simulink 6.0, and using ODE 15s variable step solver with a relative tolerance of 1eÀ3) making it difficult to integrate such a model into a larger system model.As a result, a sensitivity analysis was conducted to  Hydride Tank Aspect Ratio (height/radius) Axial Temperature Variation (K) Fig. 7 -Axial temperature variation versus tank aspect ratio.
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y 3 3 ( 2 0 0 8 ) 5 5 9 6 -5 6 0 5 identify the controlling physical features of the model and to identify techniques for reducing the model size and complexity to minimize computing resources.
The forty nodes comprise 120 integrations in the detailed model.An obvious tactic is to reduce the number of nodes in the model.Fig. 13 shows results for five cases in the sensitivity analysis regarding number of model nodes.The results from the cases having 30 and 40 nodes are nearly identical.The 30node model predicts only 0.04 g more hydrogen absorption than the 40-node model, and the temperature profiles are essentially identical.When the model size is further reduced to 20 nodes, the temperature profile remains very close to that of the 40-node model, but the hydrogen absorbed is overpredicted by a little more than 0.1 g.The temperature profile and hydrogen absorption deviate even more as the model discretization is further reduced to 10 or 5 nodes.These results suggest that the model can be reduced to 30 nodes (90 integration steps), but not further, without compromising accuracy.This reduction in size reduces the runtime to 0.07 s of real time per second modeled.
Following a suggestion of Jemni and Nasrallah [6], the analysis further investigated the sensitivity of results to the physics of heat transfer between the hydrogen gas and the solid metal hydride.Fig. 14 shows the temperature difference between hydrogen and solid metal hydride during a fill event (with hydrogen energy conservation).It is clear that the temperature difference between the phases is very small indicating a small driving force for heat transfer.If one eliminates the physics of heat transfer between the hydrogen gas and metal hydride, the dynamic energy conservation equation for the gas is simplified as follows (noting that the model now accounts for only radial transport per the previous nodal reductions): This simplification does not reduce the number of integrations in the model, but does reduce the computational time, nonetheless.To determine the model sensitivity to hydrogen energy conservation altogether one can completely eliminate hydrogen energy conservation in the dynamic model as follows: By eliminating gas energy conservation, an assumption is made that the gas in a particular node is the same temperature as the solid within that node, and that thermal energy given to or taken from the gas is negligible.This assumption reduces the number of model integration steps to 60 and produces nearly identical temperature and mass absorption levels for all conditions tested.This result is not surprising given the large difference between hydrogen heat capacity and energy of absorption.
This sensitivity analysis result can be further verified by comparing the energy that would be required to raise (or lower, upon tank emptying) the incoming gas temperature to the energy released during the exothermic reaction.The maximum internal tank temperature is determined by the particular metal hydride used and the operating pressure.For Hydralloy C5 and an inlet pressure of 13.6 atm, the van't Hoff equation (Eq.( 8)) gives the maximum temperature as 318.2 K.The change in hydrogen enthalpy from 298 K to 318.2 K can be determined by: where c p is a function of temperature [18].The value of Dh is 618 J/mol.From the previous discussion, the energy released or needed for the hydriding reaction is 28,800 J/mol.Consequently, only 2.1% of the total energy is accounted for in hydrogen temperature and can be ignored for higher level system modeling.Eliminating hydrogen energy conservation further reduces model runtime to approximately 0.005 s of real time per second of time modeled.Finally, the sensitivity analysis investigated the timescale differences between heat transfer and gas mass transfer.By inserting the values given in Table 2 into the Darcy's law equation (Eq.( 2)), one can see that hydrogen can travel the length of the 0.114 m long tank in just over 0.01 s (0.09 s/m).Heat transfer within the tank, however, takes place primarily in the radial direction.A heat transfer timescale can be obtained by dividing the mass and specific heat by the heat transfer rate, Q: Inserting property values from Table 3 gives, mc p Q ¼ ð0:435Þ ð418:7Þ ð1:0Þ ð0:0245Þ ð320 À 298Þ ¼ 338 ðs=KÞ Because the timescale of hydrogen flow is roughly four orders of magnitude faster than the heat transfer timescale, heat transfer physics dominate the dynamic response of the tank.The assumption can then be made that gas mass transfer is instantaneous, which leads to the simplification that: This in turn simplifies Eq. ( 1) to: This simplifying assumption eliminates another integration step for each node because _ m reaction is not a function of time.The sensitivity analyses suggest that several simplifying assumptions are allowable, which cumulative effect results in just one integration step per node: This further reduces the computational runtime to approximately 0.0026 s per second of real time simulated.

Results of model simplification
The hydride tank model developed herein has been incorporated into a complete RFC system model fully described   elsewhere [19].The system consists of first principle models of a 5.5-kW proton exchange membrane (PEM) fuel cell, 3.6 kW PEM electrolyzer, and three cylindrical LaNi 5 metal hydride tanks containing a total of 10,000 L of hydrogen.All three major components are linked by a common cooling system.The tank is cooled by a coolant stream and radiator when being charged and heated by the fuel cell waste heat during discharge.The ability to capture the governing physics of such a complex system in high detail with modest computing power is only possible due to the hydride tank model simplifications outlined in this work.
As an example of the effects of the hydride simplifications on the system model, Fig. 15 shows the hydrogen line pressure between the hydride tank and the fuel cell during constant fuel cell operation at 25 amp (constant fuel usage) for the complete RFC system using a 10-node hydride tank and a 30-node tank.The hydrogen line is a control volume with hydrogen from the tank as the only input and hydrogen to the fuel cell as the only output.For both cases, the tank pressure initially drops as hydrogen is desorbed, reducing tank temperature and flow rate into the line.The pressure then rises as fuel cell waste heat aids the endothermic desorption reaction resulting in more hydrogen entering the line than is used by the fuel cell.A maximum pressure is reached, after which the line pressure steadily drops as the fuel cell waste heat is no longer sufficient to release excess hydrogen from the tank due to the decreasing tank fill level.The simulation is stopped after 700 min when the tank can no longer supply hydrogen at 1 atm as required by the fuel cell.
Fig. 15 clearly shows the performance difference between the 10-node hydride tank model containing all of the physics outlined herein, and the 30-node simplified tank model.Both models contain 30 integration steps, but the 30-node model produces much smoother, more physically plausible, results.The distinct pressure fluctuations observed in the 10-node model can be traced to the temperatures of the tank nodes.Fig. 16 shows the temperature profiles for each of the 10 tank nodes during the same fuel cell operating scenario.The outermost node (Node 10) shows a smooth temperature increase due to its close proximately to the heated fuel cell waste stream.However, each subsequent node within the tank shows greater temperature fluctuations due to the thermal isolation incurred as a result of the poor thermal conductivity of the hydride alloy.
The majority of the hydrogen flowing to the fuel cell initially comes from Node 10 because the circulating heated fluid affects it most directly.As the flow rate form Node 10 wanes, Node 9 begins to empty at a much higher rate.The cooling effect associated with this produces the first temperature fluctuation observed at approximately 150 min in Fig. 16.The trend continues with each node.The 30-node model contains much smaller nodes resulting in smoother transitions, ultimately leading to the smooth gas release seen in Fig. 15.

Conclusions
A dynamic nodal model of a metal hydride hydrogen storage tank has been developed and evaluated against measured performance of a seasoned cylindrical storage tank.Comparison to data shows that the model and approach well approximate observations.As such, this proven dynamic model is shown to be useful for evaluating practical metal hydride tank design and performance.The dynamic behavior of the tank operation is then shown to be predictable with a much simpler model.The model nodal size can be reduced in combination with elimination of hydrogen energy and transport calculations to result in an accurate model containing just 30 integration steps.
The current dynamic model can assist in tank design for a particular application leading to improvements in performance.The dynamic simulation approach accounts for all of the important physics as applied to a particular geometrical configuration, yet it is simple enough to apply to integrated hydride tank -fuel cell systems in order to develop control systems and strategies.
systems.Following the original work of Mayer et al. [Mayer U, Groll M, Supper W. Heat and mass transfer in metal hydride reaction beds: experimental and theoretical results.Journal of the Less-Common Metals 1987;131:235-44], this model employs first principles heat transfer and fluid flow mechanisms together with empirically derived reaction kinetics.

ª
2008 International Association for Hydrogen Energy.Published by Elsevier Ltd.All rights reserved.

Fig. 5 -
Fig. 5 -Tank internal temperatures for each node in the radial and axial directions.

Fig. 8 -Fig. 9 -
Fig. 8 -Comparison between experimental and model temperature results for air-cooled tank using reference values for enthalpy and entropy of reaction.

Fig. 11 -
Fig. 11 -Comparison between experimental and model temperature results for air-cooled tank using experimental values for enthalpy and entropy of reaction.

Fig. 10 -Fig. 12 -Fig. 13 -
Fig. 10 -Plot showing extrapolation of experimental data used to determine enthalpy and entropy change for the current hydride reaction.

Fig. 14 -
Fig.14-Temperature difference between hydrogen gas and solid hydride when accounting for gas energy conservation.

Fig. 15 -Fig. 16 -
Fig. 15 -Plot showing line pressure between the hydride tank and a fuel cell in a complete RFC system model for two versions of hydride tank model.
Interestingly, model results showed significant error compared to the experiment when using these literature values.It is unclear if this discrepancy is due to prior poisoning of the metal alloy, or to some other degradation phenomenon.This result clearly shows that tank degradation must be considered when designing actual systems, such as RFCs, that rely on metal hydride storage.Model results were greatly improved by determining entropy and enthalpy values experimentally.

Table 1 -
Model parameter values and sources

Table 2 -
Property values used to access gas flow timescale

Table 3 -
Property values used to access heat transfer timescale