On the Kinetic Mechanisms of the Reduction and Oxidation Reactions of Iron Oxide/Iron Pellets for a Hydrogen Storage Process

Abstract

This work aims at investigating the kinetic mechanisms of the reduction/oxidation (redox) reactions of iron oxide/iron pellets under different operating conditions. The reaction principle is the basis of a thermochemical hydrogen storage system. To simulate the charging phase, a single pellet consisting of iron oxide (90% Fe₂O₃, 10% stabilising cement) is reduced with different hydrogen (H₂) concentrations at temperatures between 600 and $800 °\mathrm{C}$. The discharge phase is initiated by the oxidation of the previously reduced pellet by water vapour (H₂O) at different concentrations in the same temperature range. In both reactions, nitrogen (N₂) is used as a carrier gas. The redox reactions have been experimentally measured in a thermogravimetric analyser (TGA) at a flow rate of $250 mL/min$. An extensive literature review has been conducted on the existing reactions’ kinetic mechanisms along with their applicability to describe the obtained results. It turned out that the measured kinetic results can be excellently described with the so-called shrinking core model. Using the geometrical contracting sphere reaction mechanism model, the concentration- and temperature-dependent reduction and oxidation rates can be reproduced with a maximum deviation of less than 5%. In contrast to the reduction process, the temperature has a smaller effect on the oxidation reaction kinetics, which is attributed to 71% less activation energy ($E_{\mathrm{a},\mathrm{Re}}=56.9 \mathrm{kJ/mol}$ versus $E_{\mathrm{a},\mathrm{Ox}}=16.0 \mathrm{kJ/mol}$). The concentration of the reacting gas showed, however, an opposite trend: namely, to have an almost twofold impact on the oxidation reaction rate constant compared to the reduction rate constant.

1. Introduction

The negative effects of climate change are becoming increasingly visible in everyday life [1]. Global energy-related greenhouse and CO₂ emissions are the reasons behind this development. Despite periodic flattening, these have increased by an average of 1% per year over the last ten years [2]. To counteract this development, constructive approaches are needed to enable a global transformation of fossil fuel-dominated energy systems into energy- and resource-saving systems. A decentralised “green” energy supply system can be realised by a smart combination of wind, solar and/or bioenergy with efficient energy storage and effective energy recovery solutions. Such a sustainable energy system will enable the transition of various energy sectors such as industrial plants and urban residential areas into green societies [3,4]. The energy carrier hydrogen offers ideal conditions for such a cross-sector “green” energy supply. It can make a substantial contribution to the decarbonisation of stationary energy applications as well as the transport sector [4,5]. Additionally, H₂, as a clean energy carrier, has the potential to provide the much-needed flexibility to power systems, acting as a buffer to non-dispatchable renewable energy generation, and therefore, it can play a key role in the seasonal balancing of renewable generation by storing large amounts of it [2].

For the wide application of this potential energy carrier, suitable storage solutions are required. One storage solution is the thermochemical reduction/oxidation of metal oxides ($M_n{O}_m$)/metals (M) [6,7,8], which can be generally described as follows:

$$M_n{O}_m+m{H}_2\rightleftharpoons nM+m{H}_2\mathrm{O}.$$

(1)

While the reacting solid material does not provide a carrier material for hydrogen, it serves as a kind of catalyst for the oxidation of hydrogen to steam during the charging phase and, later on, for the reduction of steam to hydrogen during the discharging phase. Not only Fe but also other metals, such as Ce, Cu, Mn, Ni, and W, have been tested as reaction masses [9]. Based on a holistic view of the material parameters—availability, costs, reaction kinetics and energy storage density—iron proved to be the most suitable oxide carrier material [6,9,10,11].

At the beginning of the 20th century, the fundamentals of the iron redox process were used for the production of hydrogen, whereby the process was also given the name “steam–iron process” [12,13]. Nowadays, the reduction of iron oxide is mainly used in iron production [14,15]. Within the so-called direct reduction process, iron ore is reduced to metallic iron by means of a H₂/CO gas mixture. In view of the high economic importance of iron, the reduction behavior of iron oxides is a frequently studied topic in ferrous metallurgy [14,15,16].

In general, the steam iron process can be considered as a heterogeneous non-catalytic gas-solid reaction [17], in which the charging and discharging phases must be considered separately. The endothermic charging phase, in which the reduction of magnetite (Fe₃O₄) to Fe takes place, can be described as follows:

$${Fe}_3{\mathrm{O}}_4+4H_2\to 3Fe+4H_2\mathrm{O}.$$

(2)

Hydrogen reduces Fe₃O₄ first to wustite, simplified as FeO, and then to Fe. Hereby, H₂ is oxidised to steam. The corresponding oxidation reaction, which follows also the two-stage reaction principle, is given by:

$$3Fe+4H_2\mathrm{O}\to {Fe}_3{\mathrm{O}}_4+4H_2.$$

(3)

If the mass ratio of H₂ to Fe in Equations (2) and (3) is considered, a theoretical maximum storage of H₂ in Fe of $4.8 wt\%$ can be estimated, which corresponds to a gravimetric energy storage density (ESD) of $1.6 \mathrm{kwh}/(\mathrm{kg}\mathrm{Fe})$. Considering a tubular reactor being filled with spherical storage pellets each having a diameter of $d=12 \mathrm{mm}$, and made of a mixture of $90 wt\%$ hematit and $10 wt\%$ stabilisation material, the fixed-bed porosity can be estimated to 36%. The resulting volumetric ESD is then estimated to be $2.8 \mathrm{kwh}/(Lbulk)$ [18].

This exceptionally high ESD, together with the cost advantage of the storage material, offers an excellent basis for the development of an efficient hydrogen storage unit. Accordingly, several research groups have been conducting research and development tasks on this specific field. A major focus of previous research is the intermediate storage of hydrogen in a rechargeable oxide battery (ROB) with an operating temperature of $800 °\mathrm{C}$ [19,20]. By applying 30% of 8 YSZ to the iron oxide form, the highest number of 200 stable cycles, each $70 \min$ long, was obtained. However, these promising results have been realised at the cost of a remarkable reduction in the ESD from 1340 to $900 \mathrm{Wh/kg}$ due to the excessive amount of additives required [21].

Furthermore, a reformer steam iron cycle, developed by an Austrian research group, enables the production and storage of hydrogen out of bio-resources such as biogas or gasified biomass at an operating temperature of $800 °\mathrm{C}$ [9]. A $10 k{W}_{\mathrm{th}}$ laboratory prototype, which is filled in with $18 \mathrm{kg}$ of an oxygen carrier material, was presented [22]. The oxygen carrier material was applied in powder form and consists of 80% Fe₂O₃ and 20% Al₂O₃. More than 20 redox cycles could be performed in this prototype. With closing the system during the oxidation phase, the researcher showed with another test rig the possibility of producing a pressurised hydrogen up to 95 bar [23]. A further research group from the Tokyo Institute of Technology has also studied the effect of applying different additives such as Al, Cr, Cu, Ni, Rh and Zr on the storage medium in a temperature range between 300 and $380 °\mathrm{C}$ [24,25]. Iron oxide samples with a mass of $0.2 g$ (each comprising $5 \mathrm{mol\%}$ of one additive) were processed in a quartz tube reactor (length = $60 \mathrm{cm}$, i.d. = $200 \mathrm{mm}$) with reduction mixtures of CO and H₂ as well as methane. Steam with argon as a carrier gas was used as an oxidation gas. A total of 13 cycles were studied, in which Cr cations have demonstrated a prevention of the sintering of iron species during the redox reactions. In addition, Cu and Ni species as well as Rh were proven to be effective additives for hydrogen formation by the oxidation of iron with steam.

It can be concluded that in sum, the redox storage of H₂ is a promising concept due to the strong bonding of hydrogen and the high ESD [26,27]. In addition, the hydrogen is not directly stored, which is a substantial safety advantage if compared with compressed or liquified hydrogen. Furthermore, the redox storage principle offers promising integration opportunities with industrial processes operating at elevated temperatures, such as hydrogen pyrolysis [28] or the generation/combustion of H₂ using solid oxide fuel cells [29].

The overarching focus of our research group is to develop a redox storage technology based on stable iron oxide pellets [30] as a stand-alone renewable energy storage and supply system for urban quarters and industrial facilities. More about the storage process can be read elsewhere [18].

In this work, the kinetics of both reduction and oxidation (redox) reactions on an iron oxide/iron pellet have been experimentally investigated under different isothermal conditions and different concentrations of the respective reacting gas. In addition, a thorough literature review has been carried out on the existing reactions’ kinetic models and their applicability to describe the obtained experimental results. The work is devoted to gain a more systematic understanding of the kinetics of the steam-iron reactions according to (Equations (2) and (3)). To this aim, the existing kinetic models for both reduction and oxidation phases are presented and examined against our own experimental data. Both temperature and concentration dependencies of the velocity rate constants have been estimated and two conversion rate correlations for subsequent modelling activities have been suggested for both reduction and oxidation reactions on the investigated iron oxide pellets.

2. Kinetic Models

Kinetic models of chemical reactions are fundamental for the design, scale-up, optimisation, and industrial application of the reactor process. They depend on the reaction conditions such as temperature, material, particle size, flow rate, concentration, etc. For the formulation of a useful model, the prevailing chemical reaction must be correctly determined. There are many ways of classifying chemical reactions. Probably the most widely used scheme is the subdivision according to the number and type of phases involved. Generally, a differentiation is made between homogeneous and heterogeneous systems. If the reaction takes place in only one phase, it is called homogeneous. In case there are at least two phases, the reaction is termed as heterogeneous.

Additionally, the reaction can be specified according to the phase of the reactants, for example fluid–fluid or solid–fluid. On top of these classifications, the reactions can be of catalytic or non-catalytic type. Catalysts, as foreign agents that are neither reactants nor products, have a significant influence on the reaction rate. They act as intermediaries either to accelerate or to hinder the reaction. Non-catalytic reactions take place without additional substances [17,31].

In general, solid–fluid, especially solid–gas reactions are very important in many chemical and metallurgical processes, such as the reduction of oxide ores, coal gasification, roasting of sulphides and the like [14,32,33]. To study their kinetics, it is essential to define a corresponding mathematical model based on physical and chemical principles. To describe these models, a quasi-steady state is assumed, which means that the reaction is assumed to consume all the gaseous reactant supplied to the solid surface. Generally this is described with the equation:

$$\mathrm{aA}\left(g\right)+\mathrm{bB}\left(s\right)\to \mathrm{cC}\left(s\right)+\mathrm{dD}\left(g\right)$$

(4)

Here, A(g) is the reactive gas, B(s) is the active chemical solid component, C(s) is the solid product and D(g) is the gaseous product. The stoichiometric coefficients for reactants and products are given with a, b, c and d.

In order to describe the basics for the kinetic mechanisms discussed later, only the reduction of iron oxide from the reversible redox process is described in more detail below. A schematic representation of the sequential steps taking place during the reduction of a porous iron oxide particle can be seen in Figure 1. Here, the reducing gas flows through the cracks or open pores in the pellet and forms individual iron oxide regions, thus creating a flowing boundary layer. Mass transfer between the moving gas phase and the oxide surface takes place by transporting the reacting gases through this boundary layer, the so-called film layer. With the progress of the reduction phase, the reducing gas diffuses through the already reacted zone, also known as macro (1) and micro (2) pores. At the solid/gas interface, the iron oxide reduction reaction (3) described by Equation (2) takes place. This includes the adsorption of the reducing gas, the removal of oxygen from the crystal lattice (4), the formation and growth of nuclei of the reaction products wuestite or iron and finally the desorption of the oxidised gas molecules (H₂O) from the surface of the solid. The oxidised gas exits the pellet in the same way (5,6). If the kinetic rate-determining reaction is in the scenarios (1,2,5,6), it is a diffusion-limited reaction, making the overall reaction diffusion- or transport-controlled. However, if it is in the scenarios (3,4), the process becomes reaction-controlled [17,32]. A more detailed description of the physical and chemical steps taking place during such a reaction can be read in [32,34,35].

In the following subsections, the structural–physical principles utilised for the detailed description of the dynamic behavior of the reacting pellet are first derived, which is followed by the methodology of calculating the temperature- and concentration-dependent kinetic parameters.

2.1. Geometry Model

For the solid–gas system, the general molar balance equation for the reactants j = A, B (see Equation (4)) can be formulated as follows [34,36]:

$$\frac{\partial n_{\mathrm{j}}}{\partial t}=\left[\nabla \left(u{c}_{\mathrm{j}}\right)+\nabla \left(D\nabla c_{\mathrm{j}}\right)\pm R_{\mathrm{j}}\right]·\mathrm{d}{V}_{\mathrm{p}}$$

(5)

The diffusive mass transport of the reactants can be expressed by the vector of the diffusive mass flow density, which is mathematically and physically defined as the product of the molecular diffusion coefficient and the vector of the concentration gradients of the species j $\left(\nabla \left(D\nabla c_{\mathrm{j}}\right)\right)$, which is also called the 1st Fick’s law of diffusion [37]. Moreover, since the inflowing and outflowing mole streams are the same (4 moles of H₂ against 4 moles of H₂O as described by Equation (2)), the convection term $\nabla \left(u{c}_{\mathrm{j}}\right)$ can be neglected. The mass of species j formed by the chemical reaction in the balance volume is described by the rate of change of the mass $R_{\mathrm{j}}$.

Assuming a constant radius for the spherical pellet during the reaction $\left(\frac{\mathrm{d}{r}_{\mathrm{p}}}{\mathrm{d}t}=0\right)$, the differential of a porous spherical shell volume can be described as $\mathrm{d}{V}_{\mathrm{p}}={\epsilon }_{\mathrm{p}}4\pi r_{\mathrm{p}}^2\mathrm{d}{r}_{\mathrm{p}}$. Here, ${\epsilon }_{\mathrm{p}}$ specifies the volume fraction that the respective reactant j occupies on average inside the reaction volume. It is then possible to divide Equation (5) by $\mathrm{d}{V}_{\mathrm{p}}$ and add the concentration $c_{\mathrm{j}}=\frac{n_{\mathrm{j}}}{V_{\mathrm{p}}}$ to the accumulation term:

$$\frac{\partial c_{\mathrm{j}}}{\partial t}=\nabla \left(D\nabla c_{\mathrm{j}}\right)+R_{\mathrm{j}}$$

(6)

With the assumption of a constant effective diffusion coefficient $D_{\mathrm{e}}$ of the gas A through the product layer B, the equation can be reduced to:

$$\frac{\partial c_{\mathrm{A}}}{\partial t}=D_{\mathrm{e}}\left(\Delta c_{\mathrm{A}}\right)+R_{\mathrm{A}}$$

(7)

The variation of mass in the gas phase according to

$$H_2+{\mathrm{O}}^{2-}\to H_2\mathrm{O}+2e^{-}$$

(8)

is due to the mass change of the solid phase according to

$$\begin{array}{cc}\hfill \begin{array}{c}\hfill {\mathrm{Fe}}^3{\mathrm{O}}^4\end{array}+\begin{array}{c}\hfill {\mathrm{Fe}}^{2+}\end{array}+2e^{-}& \to 4\begin{array}{c}\hfill \mathrm{FeO}\end{array}\hfill \end{array}$$

(9)

$$\begin{array}{cc}\hfill \begin{array}{c}\hfill \mathrm{FeO}\end{array}& \to \begin{array}{c}\hfill \mathrm{Fe}\end{array}+\begin{array}{c}\hfill \mathrm{O}\end{array}.\hfill \end{array}$$

(10)

Including the stochiometric coefficients according to Equation (4), the mass change of the solid B can thus be described, assuming irreversibility and first-order rate in the concentration of the gas reactant, as follows [38,39]:

$$\frac{\mathrm{d}{c}_{\mathrm{B}}}{\mathrm{d}t}=-b{R}_{\mathrm{B}}$$

(11)

with

$$R_{\mathrm{B}}=f\left(\alpha \right)·c_{\mathrm{A}},$$

(12)

where $c_{\mathrm{A}}$ is the reaction gas concentration and $f\left(\alpha \right)$ describes the amount of solid B per unit volume of the solid. This function can be expressed after meeting some assumptions concerning the behavior of solid matrix during the conducted reaction. In the literature [40,41,42,43], the following four different functional expressions are present:

• In Reaction order models, the reaction rate is proportional to the concentration, where the remaining amount or fraction of the reactant is raised to a particular power (integral or fractional), which represents the reaction order $n_{\mathrm{RB}}$. The general formulation for a reaction-based model is:

  $$f\left(\alpha \right)={(1-\alpha )}^{n_{\mathrm{RB}}}$$

  (13)
• Geometric contraction models assume that the adsorption and desorption of the reaction gas occurs uniformly on the surface of the particle. The reaction rate is controlled by the resulting expansion of the contact area towards the centre of the reactant. Assuming that all particles have the same shape and size, this model is easy to parameterise and shows good agreement for metal oxide/metal redox systems. The two relevant kinetic parameters, the generalised rate constant and the morphology $n_{\mathrm{GC}}$ (cubic, cylindrical, spherical), can be determined by fitting a single isothermal curve. In general, this can be described with the following algebraic equation:

  $$f\left(\alpha \right)=(n_{\mathrm{GC}}+1){(1-\alpha )}^{n_{\mathrm{GC}}/(n_{\mathrm{GC}}+1)}$$

  (14)
• Diffusion models describe complex concentration equilibrium processes through porous solid matrices or simply pore diffusion, which can play an important role in gas–solid reactions. If the reacting solid B has a porosity ($1-\epsilon$), the description of diffusion through the micro- and macro-pore volume is necessary for the reacting gas A to access the surface of the solid B; the removal of the gaseous product D also occurs in the same way. However, pore diffusion can also be an important component in the reaction of non-porous solids. If the solid product layer C formed is itself porous, the supply of the gaseous reactant A and the removal of the gaseous product D can take place by diffusion through the porous product layer. Depending on the model used, e.g., Gingstling–Brounshtein, Jander or Chou, the mathematical description differs and will be considered in more detail later.
• Nucleation and nuclei growth models describe the formation and growth of nuclei, which are finite quantities of the product inside the reactant lattice. This includes crystallisation, crystallographic transition, decomposition, adsorption, hydration and desolvation. In general, the rate-controlling step can be the solid diffusion or movement of the interface.

For modelling the structural change of the reaction mass inside the pellet depending on the reaction progress of the heterogeneous, non-catalytic redox reaction, various models can be found in the literature [17,37,44]. Among these models, it can be distinguished between porous and non-porous reaction solid models. Depending on the structural model used, some of these functional mechanisms are applied or even combined, which implies the selection of the appropriate initial and boundary conditions for Equation (7). To formulate an adequate model for describing the complex structure and, consequently, the behavior of the iron oxide pellet during the reaction, the following assumptions deem necessary:

• Isothermal condition dominates $\mathrm{d}T/\mathrm{d}t=0$;
• Based on Equations (2) and (3), equimolar counterdiffusion is valid $n_{\mathrm{A}}=-n_{\mathrm{D}}$;
• Reactions proceed independently and without interaction between the gaseous species;
• No crack formation occurs during the reaction;
• Pressure is uniform inside and around the pellets;
• Pellet has a uniform and constant porosity $\mathrm{d}\epsilon /\mathrm{d}t=0$;
• Chemical reaction-controlled process inside each grain is reversible and proceeds topochemically.

The simplest structural model, described in this work, is the progressive-conversion model (PCM). PCM assumes that, according to Equation (4), the reaction gas (A) enters the solid (B) at any time and reacts there. If the diffusion rate through the pores of the particle is much higher than the reaction rate $D_{\mathrm{e}}>>R_{\mathrm{A}}$, the concentration profile of the reaction gas in Equation (7) will be uniformly distributed over the particle $\frac{{\partial }^2{c}_{\mathrm{A}}}{\partial x_i^2}=0$.

Consequently, the chemical reaction $R_{\mathrm{A}}$ is the only step in the PCM that controls the molar mass change, as illustrated in Figure 2a. Assuming that the concentration of the reaction gas is so high that it remains unchanged during the process $\frac{\mathrm{d}{c}_{\mathrm{A}}}{\mathrm{d}t}=0$, and assuming a first-order reaction takes place inside the solid B, the molar mass consumption of component B can be described (see Equation (11)) as [17]:

$$\frac{\mathrm{d}{n}_{\mathrm{B}}}{\mathrm{d}t}=b\mathrm{d}{V}_{\mathrm{p}}k{c}_{\mathrm{A}}=\epsilon 4\pi r_{\mathrm{p}}^{2}bk{c}_{\mathrm{A}}\mathrm{d}{r}_{\mathrm{p}}$$

(15)

The analytical solution for the consumption of component B (${\alpha }_{\mathrm{B},\mathrm{PCM}}\left(t\right)$) of this equation is:

$${\alpha }_{\mathrm{B},\mathrm{PCM}}=\frac{t}{{\tau }_{\mathrm{PCM}}}\wedge {\tau }_{\mathrm{PCM}}=\frac{{\rho }_{\mathrm{B}}}{b{a}_{\mathrm{V}}{k}_{\mathrm{PCM}}{c}_{\mathrm{A}}}.$$

(16)

where ${\rho }_{\mathrm{B}}$ is the molar density of the reactive component in the solid B, $k_{\mathrm{PCM}}$ is the reaction rate constant, $c_{\mathrm{A}}$ is the molar concentration of the reaction gas, $a_{\mathrm{V}}$ is the specific surface of the particle (area per unit volume) which remains constant during the proces and ${\tau }_{\mathrm{PCM}}$ is the time needed for completing the reaction. Assuming a constant reaction gas concentration inside the entire particle (equal to that in the bulk) and a first-order reaction at the surface of the reaction material, the PCM is the choice to describe leaching processes [44].

In order, however, to model the kinetics of the redox reactions, most of the literature reported that a concentration change over the pellet radius must be assumed $\frac{\mathrm{d}{c}_{\mathrm{A}}\left(t\right)}{\mathrm{d}r}\ne 0$ [17,33,34]. Having a metal oxide reacting with a reactant gas, the reaction zone changes or moves from the outer shell of the solid material towards its interior sections during the reaction and, consequently, the gas concentration does change as it continues to react with the shrinking product layer. Such a behavior can be described with the so-called shrinking-core model (SCM). SCM describes an active pellet with a changing reaction zone $A_{\mathrm{c}}=4\pi {\left(r_{\mathrm{c}}\left(t\right)\right)}^2$. The reaction takes place, first, on the outer skin $A_{\mathrm{p}}=4\pi r_{\mathrm{p}}^2$ of the particle and then moves into the solid. This leaves behind a completely transformed material as well as the inert solid material. This is usually referred to as “ash”. At any time, there is an unreacted core of material that shrinks during the reaction, as illustrated in Figure 2b [17].

The authors of [16,17], have sliced and examined the cross-section of partly reacted solid particles and reported that the unreacted solid material inside the grains are surrounded by a layer of ash. The unreacted core may not always be as sharply defined as the model assumes it; nevertheless, the literature do recommend the SCM for describing the redox behavior of metal oxides [17,33,45,46,47]. The main advantages of the SCM are the simple parameterisation and the existence of an analytical solution for the particle conversion, provided that the gas concentration surrounding the particle (inside the edges or the bulk phase) is constant, i.e., for $r>r_{\mathrm{p}}\Rightarrow \frac{\mathrm{d}{c}_{\mathrm{A}}}{\mathrm{d}t}=c_{\mathrm{A},\mathrm{Bulk}}=c_{\mathrm{A}}$. It also provides an analytical function for the global reaction rate of the entire pellet as a function of the reaction gas concentration as well as the reaction temperature. This is quite useful, especially when the particle model needs to be coupled with a reactor model that is more complex than a simple plug flow [39].

To describe the particle model, the different kinetic mechanisms must be considered. SCM consists of three kinetic control mechanisms, which results from Equation (7) [36] and are described briefly in the following paragraphs:

(1) Diffusion of the reacting gas through the gas film layer (boundary condition) at the corresponding mass transfer coefficient $\beta \left(\frac{\mathrm{m}}{\mathrm{s}}\right)$, also termed as the film theory,

$$D_{\mathrm{e}}\frac{\mathrm{d}{c}_{\mathrm{A}}\left(t\right)}{\mathrm{d}r}{|}_{r=r_{\mathrm{p}}}=\beta \left(c_{\mathrm{A}}-c_{\mathrm{A}}\left(t\right){|}_{r=r_{\mathrm{p}}}\right)t>0$$

(17)

$$\Rightarrow -\frac{\mathrm{d}{n}_{\mathrm{A}}}{\mathrm{d}t}=4\pi {\left(r_{\mathrm{c}}\left(t\right)\right)}^2\beta \left(c_{\mathrm{A}}-c_{\mathrm{A}}{|}_{r=r_{\mathrm{p}}}\right)$$

(18)

(2) Diffusion through the porous product layer (ash layer) with the effective diffusion coefficient $D_{\mathrm{e}}\left(\frac{{\mathrm{m}}^2}{\mathrm{s}}\right)$

$$-\frac{\mathrm{d}{n}_{\mathrm{A}}}{\mathrm{d}t}=4\pi {\left(r_{\mathrm{c}}\left(t\right)\right)}^2{D}_{\mathrm{e}}\frac{\mathrm{d}{c}_{\mathrm{A}}}{\mathrm{d}r}{|}_{r=r_{\mathrm{c}}\left(t\right)}$$

(19)

(3) Chemical reaction for an irreversible first-order reaction with the kinetic rate constant $k\left(\frac{\mathrm{m}}{\mathrm{s}}\right)$

$$-\frac{\mathrm{d}{n}_{\mathrm{A}}}{\mathrm{d}t}=4\pi {\left(r_{\mathrm{c}}\left(t\right)\right)}^{2}k{c}_{\mathrm{A}}{|}_{r=r_{\mathrm{c}}\left(t\right)}.$$

(20)

The concentration gradient $\frac{\mathrm{d}{c}_{\mathrm{A}}}{\mathrm{d}r}$ in the formed product layer can be derived from the molar balance in spherical coordinates by neglecting the convection and reaction terms and assuming pseudo-stationary conditions, as follows:

$$0=\nabla c_{\mathrm{A}}=\frac{{\mathrm{d}}^2{c}_{\mathrm{A}}}{\mathrm{d}{r}^2}+\frac{2}{r}\frac{\mathrm{d}{c}_{\mathrm{A}}}{\mathrm{d}r}$$

(21)

Equation (21) is also known as Fick’s second law of diffusion. Applying the boundary conditions for $r=r_{\mathrm{p}}\Rightarrow c_{\mathrm{A}}=c_{\mathrm{A}}{|}_{r=r_{\mathrm{p}}}$ and for $r=r_{\mathrm{c}}\left(t\right)\Rightarrow c_{\mathrm{A}}=c_{\mathrm{A}}{|}_{r=r_{\mathrm{c}}\left(t\right)}$ and carrying out a twofold integration of Equation (21) results in the following expression for the radius dependency of the gas concentration:

$$c_{\mathrm{A}}\left(r\right)=c_{\mathrm{A}}{|}_{r=r_{\mathrm{c}}\left(t\right)}+\left(c_{\mathrm{A}}{|}_{r=r_{\mathrm{p}}}-c_{\mathrm{A}}{|}_{r=r_{\mathrm{c}}\left(t\right)}\right)\frac{1-\frac{r_{\mathrm{c}}\left(t\right)}{r}}{1-\frac{r_{\mathrm{c}}\left(t\right)}{r_{\mathrm{p}}}}$$

(22)

The concentration gradient at the point $r=r_{\mathrm{c}}\left(t\right)$ can be expressed as follows:

$$\frac{\mathrm{d}{c}_{\mathrm{A}}}{\mathrm{d}r}{|}_{r=r_{\mathrm{c}}\left(t\right)}=\frac{c_{\mathrm{A}}{|}_{r=r_{\mathrm{p}}}-c_{\mathrm{A}}{|}_{r=r_{\mathrm{c}}}}{r_{\mathrm{c}}\left(t\right)\left(1-\frac{r_{\mathrm{c}}\left(t\right)}{r_{\mathrm{p}}}\right)}$$

(23)

By substituting Equations (18)–(20) and reformulating it, the surface concentration $c_{\mathrm{A}}{|}_{r=r_{\mathrm{c}}\left(t\right)}$ as a function of $r_{\mathrm{c}}\left(t\right)$ is obtained as follows

$$c_{\mathrm{A}}{|}_{r=r_{\mathrm{c}}\left(t\right)}=\frac{c_{\mathrm{A}}}{1+\frac{k}{\beta }{\left(\frac{r_{\mathrm{c}}\left(t\right)}{r_{\mathrm{p}}}\right)}^2+\frac{k{r}_{\mathrm{c}}\left(t\right)}{D_{\mathrm{e}}}\left(1-\frac{r_{\mathrm{c}}\left(t\right)}{r_{\mathrm{p}}}\right)}$$

(24)

Together with Equation (20), the time rate of the converted molar mass of A for a spherical pellet with film and ash diffusion as well as chemical reaction is obtained:

$$\frac{\mathrm{d}{n}_{\mathrm{A}}}{\mathrm{d}t}=-\frac{4\pi {\left(r_{\mathrm{c}}\left(t\right)\right)}^{2}k{c}_{\mathrm{A}}}{1+\frac{k}{\beta }{\left(\frac{r_{\mathrm{c}}\left(t\right)}{r_{\mathrm{p}}}\right)}^2+\frac{k{r}_{\mathrm{c}}\left(t\right)}{D_{\mathrm{e}}}\left(1-\frac{r_{\mathrm{c}}\left(t\right)}{r_{\mathrm{p}}}\right)}$$

(25)

A structural study of the pellet requires an expression for $r_{\mathrm{c}}\left(t\right)$ as a function of time. The following correlation exists between the time derivative of the converted solid B and $r_{\mathrm{c}}\left(t\right)$:

$$\frac{\mathrm{d}{n}_{\mathrm{B}}}{\mathrm{d}t}=-\frac{4\pi {\left(r_{\mathrm{c}}\left(t\right)\right)}^2{\rho }_{\mathrm{B}}}{M_{\mathrm{B}}}\frac{\mathrm{d}{r}_{\mathrm{c}}\left(t\right)}{\mathrm{d}t}=\frac{b\mathrm{d}{n}_{\mathrm{A}}}{\mathrm{d}t}$$

(26)

Finally, the shrinking unreacted core can be described by applying the stoichiometric coefficients from Equation (4) as follows:

$$\frac{\mathrm{d}{r}_{\mathrm{c}}\left(t\right)}{\mathrm{d}t}=-\frac{c_{\mathrm{A}}}{c_{\mathrm{B}},0}\frac{bk}{1+\frac{k}{\beta }{\left(\frac{r_{\mathrm{c}}\left(t\right)}{r_{\mathrm{p}}}\right)}^2+\frac{k{r}_{\mathrm{c}}}{D_{\mathrm{e}}}\left(1-\frac{r_{\mathrm{c}}\left(t\right)}{r_{\mathrm{p}}}\right)}$$

(27)

$M_{\mathrm{B}}\left(\frac{\mathrm{k}\mathrm{g}}{\mathrm{mol}}\right)$ stays for the molar mass of solid B, $b(-)$ for its stoichiometric coefficient and ${\rho }_{\mathrm{B}}\left(\frac{\mathrm{k}\mathrm{g}}{{\mathrm{m}}^3}\right)$ for its density. If the gas inflows continuously, $c_{\mathrm{A}}$ can be assumed to be constant in time and position. With the relationship between the degree of conversion of the solid ${\alpha }_{\mathrm{B}}$ and $r_{\mathrm{c}}$

$${\alpha }_{\mathrm{B}}\left(t\right)=1-\frac{c_{\mathrm{B}}\left(t\right)}{c_{\mathrm{B}},0}=1-{\left(\frac{r_{\mathrm{c}}\left(t\right)}{r_{\mathrm{p}}}\right)}^3$$

(28)

$$\Rightarrow \frac{\mathrm{d}{\alpha }_{\mathrm{B}}}{\mathrm{d}t}=-\frac{3}{r_{\mathrm{p}}}{\left(\frac{r_{\mathrm{c}}\left(t\right)}{r_{\mathrm{p}}}\right)}^2\frac{\mathrm{d}{r}_{\mathrm{c}}\left(t\right)}{\mathrm{d}t}=\frac{3}{r_{\mathrm{p}}}\frac{1}{\frac{1}{\beta }+\frac{1}{k}{\left(1-{\alpha }_{\mathrm{B}}\right)}^{-2/3}+\frac{r_{\mathrm{p}}}{D_{\mathrm{e}}}\left[{\left(1-{\alpha }_{\mathrm{B}}\right)}^{-1/3}-1\right]}\frac{b{c}_{\mathrm{A}}}{c_{\mathrm{B}},0}$$

(29)

the function $f\left({\alpha }_{\mathrm{B}}\right)$ for the SCM is defined as follows [34,36]:

$$f\left({\alpha }_{\mathrm{B}}\right)=-\frac{1}{b}\frac{c_{\mathrm{B}},0}{c_{\mathrm{A}}}\frac{\mathrm{d}{\alpha }_{\mathrm{B}}}{\mathrm{d}t}$$

(30)

As long as an adequate flow around the pellet is provided, the mass transfer $\beta >>k$. Applying this condition in Equation (31) results in:

$$\Rightarrow f\left({\alpha }_{\mathrm{B}}\right)=\frac{-3}{\frac{r_{\mathrm{p}}}{k}{\left(1-{\alpha }_{\mathrm{B}}\right)}^{-2/3}+\frac{r_{\mathrm{p}}^2}{D_{\mathrm{e}}}\left[{\left(1-{\alpha }_{\mathrm{B}}\right)}^{-1/3}-1\right]}$$

(31)

The overall process rate is determined by the slowest of the two sequential mechanisms, which take place one after the other during the reaction. In other words, if the diffusion through the porous solid layer is the rate-controlling mechanism, the following applies:

$$\frac{1}{k}\to 0\Rightarrow f\left({\alpha }_{\mathrm{B}}\right)=\frac{-D_{\mathrm{e}}}{r_{\mathrm{p}}^2}\frac{-3}{{\left(1-{\alpha }_{\mathrm{B}}\right)}^{-1/3}-1}$$

(32)

In case the chemical reaction is considered to be the rate-controlling mechanism

$$\frac{1}{D_{\mathrm{e}}}\to 0\Rightarrow f\left({\alpha }_{\mathrm{B}}\right)=\frac{-k}{r_{\mathrm{p}}}3{\left(1-{\alpha }_{\mathrm{B}}\right)}^{2/3}$$

(33)

With the effective diffusion coefficient $D_{\mathrm{e}}$, it is possible to describe the resistance of the ash macro pores to the fluid flow in a simplified way. In general, two types of the pores must be considered: pores that are formed during the reaction and those that are available in the raw material.

Different sizes and shapes of these pores lead to different diffusion mechanisms for the reaction gas flow [48]. If the mass transport in the marco pores is the rate-limiting step, it is possible to describe the reaction process using the SCM. In case of high porous masses, which consists of small particles, the effect of gas diffusion as well as the structural properties, such as diameter, porosity and tortuosity, becomes more important. The sharp boundary between the reacted and the unreacted zone can no longer be confirmed for such reaction masses [49,50]. In order to integrate additional reaction zones of particles in the porous mass, the well-known grain model (GM) is applied.

With the GM, not only the diffusion through the formed product layer is taken into account but also the structure of the unreacted solid particles. The gas penetrates through the pores and reacts on the grains’ surface, as shown in Figure 1 and Figure 2c. According to Szekely [51], the solid reactant comprises micro spherical grains, with a reacting core that shrinks as long as the reaction continues. These grains are assumed to be dense and fine with a constant radius $r_{\mathrm{g}}$. Structures of this type can be observed in many porous materials [35,52]. Complementary to the derivation of the SCM in Equation (31), which describes the reaction behavior of one grain, a finite number of micro grains is integrated into the equations to describe the GM, which increases, however, the numerical effort. The reader is referred to [34,39,49,53] for a detailed derivation of the GM.

So far, a general description of the reaction behavior on one pellet was derived, on which basis the prevailing dominant mechanism can be determined, i.e., whether the process is diffusion or chemical reaction dominated. In each case, the relevant kinetic parameters can be estimated.

2.2. Estimation of the Kinetic Parameters

In order to investigate the influences of the different reaction conditions, a number of independent algebraic equations have been presented in the literature, which are useful for describing the process within a certain range of reaction conditions [40,54,55,56,57]. The basis for this is the general description of the solid-state reaction rate:

$$\frac{d\alpha }{dt}=k_{\mathrm{A}}·f\left(\alpha \right),$$

(34)

where the concentration and temperature dependence is described by the modified Arrhenius rate constant $k_{\mathrm{A}}$$\left(\frac{1}{s}\right)$, and the degree of conversion is defined $f\left(\alpha \right)(-)$ by an algebraic equation, which represents a physical approximation of the kinetic model and is already discussed in Section 2.1. Exchanging $k_{\mathrm{A}}$ by the modified Arrhenius equation results in:

$$\frac{d\alpha }{dt}=\underset{k_{\mathrm{A}}}{\underbrace{A{e}^{-\frac{E_{\mathrm{a}}}{RT}}{C}^n}}·f\left(\alpha \right)$$

(35)

with $\alpha$ as the reacted fraction at time t, $E_{\mathrm{a}}$ as the activation energy, A as a pre-exponential factor of the Arrhenius relationship, and T as the absolute temperature. In addition, the reaction gas concentration C appears, which is raised to an empirical exponent n [58,59,60], which has to be determined by fitting the experimental data to Equation (35).

The algebraic equation $f\left(\alpha \right)$ is generally classified based on the graphical shape of their isothermal curves or on their mechanistic assumptions. Based on the shape of the isothermal curves of ($d\alpha /dt$ vs. $\alpha$), the kinetic models can be grouped into acceleratory, deceleratory, constant/linear, or sigmoidal models (see Figure 3). Acceleratory models are those in which the reaction rate ($d\alpha /dt$) is increasing as the reaction proceeds (nucleation models). Similarly, deceleratory reaction rates decrease with the reaction progress (phase boundary controlled reaction models or diffusion models). The reaction rate remains constant for linear models (reaction order-based models), while sigmoidal models show a bell-shaped relationship (nuclei growth models).

The mechanistic assumptions $f\left(\alpha \right)$ evaluated from the plot $\alpha$ vs. t are divided into the already discussed physical processes: Nucleation and nuclei growth (mainly represented by the Avrami–Erofeev equation), progressive conversion (reaction order-based models) as well as diffusion and geometrical contraction [40,61]. The methods for graphical presentation are easy means of visually determining the most appropriate model for a particular dataset. An extraction of the most commonly used kinetic models is presented in

Table 1. Algebraic expressions for $f\left(\alpha \right)$ and $g\left(\alpha \right)$ functions for different kinetic models [40,50,62].

Reaction Mechanism	Symbol	Kinetic Mechanism Function
		Differential Form	Integral Form
		$\mathit{f}\left(\mathit{\alpha }\right)=\frac{\mathbf{1}}{{\mathit{k}}_{\mathbf{A}}}\frac{\mathit{d}\mathit{\alpha }}{\mathit{dt}}$	$\mathit{g}\left(\mathit{\alpha }\right)={\mathit{k}}_{\mathbf{A}}\mathit{t}$
Nucleation models (Avarami–Erofe’ev)
Nucleation and growth, n = order of reaction	An	$n(1-\alpha ){[-ln(1-\alpha )]}^{(n-1)/n}$	${[-ln(1-\alpha )]}^{1/n}$
Phase boundary-controlled reaction
Geometrical contracting infinite slab (Film Diffusion control)	R1	1	$\alpha$
Geometrical contracting cylinder (Chemical reaction control)	R2	$2{(1-\alpha )}^{1/2}$	$1-{(1-\alpha )}^{1/2}$
Geometrical contracting sphere (Chemical reaction control)	R3	$3{(1-\alpha )}^{2/3}$	$1-{(1-\alpha )}^{1/3}$
Diffusion models
3D Diffusion Jander eqn.	D3	$\frac{3{(1-\alpha )}^{2/3}}{2(1-{(1-\alpha )}^{1/3})}$	${[1-{(1-\alpha )}^{1/3}]}^2$
Ginstling–Brounshtein (Ash Layer, Diffusion control)	D4	$\frac{3}{2\left({(1-\alpha )}^{-1/3}-1\right)}$	$1-\frac{2}{3}\alpha -{(1-\alpha )}^{2/3}$
Reaction order-based models
First-order	F1	$(1-\alpha )$	$-ln(1-\alpha )$
Second-order	F2	${(1-\alpha )}^2$	$[1/(1-\alpha \left)\right]-1$

. A more detailed description can be read in [40].

For the parametrisation of the mechanistic models, it is common to use the “model-fitting” method [62,63]. Based on the fitting method, $f\left(\alpha \right)$ as well as $k_{\mathrm{A}}$ can be determined by fitting the integral form of Equation (34), which is described by

$$\begin{array}{c}{\int }_0^{\alpha }\frac{d\alpha }{f\left(\alpha \right)}=k_{\mathrm{A}}{\int }_0^{t}dt\end{array}$$

(36)

$$\begin{array}{c}\to g\left(\alpha \right)=k_{\mathrm{A}}t\end{array}$$

(37)

against the experimental data, which represents the mass change of the sample. The integral form $g\left(\alpha \right)$ can be determined by integrating the kinetic models of $f\left(\alpha \right)$ listed in Table 1. By means of a mathematical optimisation method, the activation energy $E_{\mathrm{a}}$ can, subsequently, be determined by using the calculated rate constant k from the isothermal experimental data [57].

For the purpose of completeness, the model of Sestak and Berggren [64] should be mentioned here, which showed that it is possible to express the integral form of $g\left(\alpha \right)$ in a general analytic form that represents all models in a single general expression:

$$g\left(\alpha \right)={\alpha }^m{(1-\alpha )}^o{[-ln(1-\alpha )]}^p$$

(38)

where m, p and o are constants with the listed possible combinations in

Table 2. Combination of o, m and p as well as a description to possible mechanistic models [64]; models marked in italics are not considered here.

o	m	p	Is Suitable for:
−	−	−	Phase boundary controlled reaction
✓	−	−	Phase boundary-controlled reaction (R2, R3), Reaction order models (F1)
−	✓	−	Mampel Equation (P2, P3), linear nucleation and diffusion models
−	−	✓	Gingstling–Brounshtein (D4), 2D diffusion models
✓	✓	−	Nucleation, latter stages of linear growth of nuclei
✓	−	✓	Growth of nulcei (A1.5–A4)
−	✓	✓	Unjustified as yet
✓	✓	✓	Any complicated case, unjustified as yet

along with the most suitable application fields of each mechanistic model.

By assigning empirical values for these three variables, the listed physical models in Table 1 can be obtained [64]. An overview of more parameter values and their combinations for more models can be read elsewhere [65].

In order to evaluate the reaction models based on the measured isothermal gravimetric kinetic data, the instantaneous reduction/oxidation fraction ${\alpha }_{\mathrm{E}}\left(t\right)$ is defined by Equation (39) as the instantaneous mass change of the sample $m\left(t\right)$ relative to its initial mass $m_0$ divided by the final mass change $m_{\infty }$ reached after a given maximum reaction time.

$${\alpha }_{\mathrm{E}}\left(t\right)=\frac{m_0-m\left(t\right)}{m_0-m_{\infty }}$$

(39)

The smallest representative investigation unit of a gas–solid reaction system is the interaction of a single particle with a moving gas stream. Since the investigation of single-particle systems is more transparent and the results can be generalised, at least in principle [34], to more complex multi-particle arrangements, a single-pellet measuring device will be described in the following section.

3. Materials and Methods

For investigating the solid–gas redox reaction kinetics in this work, the so-called “engineering approach” is followed, in which a relatively big sample mass is investigated under well-defined external boundary conditions; namely, temperature and gas concentration, in order to obtain the necessary information for the later practical application of the process [50]. A commonly used methodology for such an engineering approach is the determination of the reaction kinetics by means of a thermogravimetric analyser. The set-up used to investigate the iron oxide samples is illustrated schematically in Figure 4, where a hydrogen-compliant thermogravimetric analyser (STA 449-F3, NETZSCH Gerätebau GmbH, Selb, Germany) forms the core component.

The reduction of the iron oxide pellet is triggered by flowing a mixture of nitrogen and hydrogen into the TGA device, where the sample is placed. Different gas mixtures are adjusted by means of a mass flow controller for each gas (EL-FLOW, Bronkhorst, Kamen, Germany) and fed into the TGA device. The oxidation of the previously reduced samples takes place in an H${}_2$O-H${}_2$ gas mixture atmosphere. The required steam is produced in a steam generator (DV2MK, aDROP Feuchtemeßtechnik GmbH, Fürth, Germany) and transferred to the TGA by means of a heated metallic hose to avoid condensation. Furthermore, a vacuum pump is used to evacuate the experimental setup after each experiment for safety reasons. The measured values of the experimental setup are stored locally on the measurement computer with a sampling rate of $200\mathrm{points}/ \min$. The accuracies of the applied devices are listed in

Table 3. List of all measuring device accuracies for the experimental setup.

Device	Physical Unit	Measurement Accuracy
STA 449-F3, Netzsch	Mass	35,000 mg ± 0.0001 mg
	Temperature	$1600 °\mathrm{C}\pm 0.3 °\mathrm{C}$
EL-FLOW, Bronkhorst	Gas flow (H2, N2)	$250 \mathrm{mL/min}$ ± 0.5%
Steam generator DV2MK, ADROP	Mass flow (H2O)	$25 \mathrm{g/h}$ ± 1%

.

The investigated sample, based on the recipe of [30], was acquired from the company, Nacompex GmbH, Bad Saarow, Germany. The almost spherical sample consists of a mixture of 90% hematite (Fe₂O₃) and 10% support material of high-alumina cement [66]. The main function of the support material is to ensure a high dispersion and a high mechanical strength of the pellet. After mixing and shaping processes, the pellets were dried for 28 days under ambient conditions before they were sintered at $950 °\mathrm{C}$ for $3 h$. The investigated sample has a diameter of $d_{\mathrm{s}}\approx 12.7 \mathrm{mm}$ with a total mass of $m_{\mathrm{s}}=3658.5 \mathrm{mg}$.

Before the very first reduction reaction, the raw sample consists of 90% Fe₂O₃, which is reduced to Fe using H₂. In the subsequent reactions, and because of the adjusted operating conditions (

Table 4. Influx molar concentrations of the reducing and oxidation gas mixtures during the four redox experiments ro${}_1$ to ro${}_3$.

Label	Reduction		Oxidation
	H2/N2		H2O/N2
	$\frac{\mathbf{mL}}{\mathbf{min}}$/$\frac{\mathbf{mL}}{\mathbf{min}}$	mol%/mol%	$\frac{\mathit{g}}{\mathit{h}}$/$\frac{\mathbf{mL}}{\mathbf{min}}$	mol%/mol%
ro${}_1$	63/187	25/75	1.8/162	16/84
ro${}_2$	125/125	50/50	3.6/125	33/67
ro${}_3$	245/5	98/2	5.8/60	67/33

), Fe can only oxidise to Fe₃O₄ as can be observed in Section 4.1. From there on, Fe₃O₄ is converted into Fe and vice versa. More details on this conversion can be read in [67].

As described before, during the experimental campaign for the determination of the reaction kinetics (see Figure 4), the molar concentrations for reduction and oxidation were changed according to Table 4 for the different redox runs (ro${}_1$–ro${}_3$).

The oxidation atmosphere of ro${}_3$ results from the maximum allowed steam concentration for the safe and accurate operation the STA device. At higher concentrations, undesired condensation may take place inside the weighing section of the STA.

Previous studies have reported that the iron composite material can sinter if repeated reactions are performed above an oxidation temperature of $800 °\mathrm{C}$ [18,67]. For this reason, the reduction/oxidation investigations have been carried out under the isothermal conditions of $600 °\mathrm{C}$, $700 °\mathrm{C}$, and $800 °\mathrm{C}$.

4. Results and Discussion

This section is divided into two main subsections: namely, Reduction and Oxidation. At the beginning of each subsection, an algebraic equation from Table 1 is fitted to the experimental values using the model-fitting method for isothermal experiments [41]. As previously explained, the kinetic parameters, activation energy $E_{\mathrm{a}}$ and Arrhenius’ pre-exponential factor A, are determined for the respective redox reaction by means of a mathematical optimisation fitting to the algebraic equations. To solve these non-linear curve-fitting problems (data fitting), the least squares method is applied in MATLAB R2019b [68]. Generally formulated, the optimisation problem for the velocity rate constant is

$$\underset{k_{\mathrm{A}}}{\min }{\left|\left|F_{\exp }-\frac{d\alpha \left(k_{\mathrm{A}}\right)}{dt}\right|\right|}_2^2=\underset{k_{\mathrm{A}}}{\min }\sum _i{\left(F_{\exp ,i}-\frac{d{\alpha }_i\left(k_{\mathrm{A}}\right)}{dt}\right)}^2$$

(40)

and for the activation energy is

$$\underset{A,E_{\mathrm{a}}}{\min }{\left|\left|F_{k_{\mathrm{A}},\exp }-f(A,E_{\mathrm{a}})\right|\right|}_2^2=\underset{A,E_{\mathrm{a}}}{\min }\sum _i{\left(F_{k_{\mathrm{A}},\exp ,i}-f_i(A,E_{\mathrm{a}})\right)}^2$$

(41)

with $F_{\exp }$ representing the experimental data, $\frac{d\alpha }{dt}$ as well as the searched coefficient $k_{\mathrm{A}}$ from Equation (34). $F_{k_{\mathrm{A}},\exp }$ is a function of the determined velocity rates from the experimental data and $f(A,E_{\mathrm{a}})$ with the searched parameters $A,E_{\mathrm{a}}$ for the Arrhenius formalism of Equation (35). A flowchart describing the estimation methodology for the kinetic parameters is depicted in Figure 5.

The absolute deviation (AD) between the experimentally measured reduction or oxidation fraction ${\alpha }_{\mathrm{Exp}}$ and the mathematically fitted values ${\alpha }_{\mathrm{Fit}}$ is defined by Equation (42). For better reading and interpreting the diagrams, AD is expressed in (%) rather than in a dimensionless fraction.

$$\mathrm{AD}=\left|{\alpha }_{\mathrm{Exp}}-{\alpha }_{\mathrm{Fit}}\right|·100\%$$

(42)

4.1. Reduction

As described previously, the reduction behavior of a metal oxide is a complex thermochemical gas–solid process, which has been in depth investigated with the main target to improve the theoretical models and the applied algebraic rate equations [32,34,40,69]. In this subsection, the focus is specifically on the reduction of iron oxide using hydrogen, for which a summary of the relevant literature can be found in

Table 5. Summary of researcher on the iron oxides reduction by hydrogen as well as oxidation by steam and the used kinetics models.

Reaction	Temperature	Conc.	Geometry	Model	${\mathit{E}}_{\mathbf{a}}$	Experiment	Ref.
$\begin{array}{c}\hfill {\mathrm{Fe}}_3{\mathrm{O}}_4\end{array}\to \begin{array}{c}\hfill \mathrm{Fe}\end{array}$	$1173 K$	$\begin{array}{c}\hfill H_2\end{array}/\begin{array}{c}\hfill C\mathrm{O}\end{array}=1.5$	Pellet, $\mathrm{d}=15.2 \mathrm{mm}$	R3	$75.4 \mathrm{kJ/mol}$	Isothermal	[53,72]
$\begin{array}{c}\hfill \mathrm{FeO}\end{array}\to \begin{array}{c}\hfill \mathrm{Fe}\end{array}$	$1173 K$	$\begin{array}{c}\hfill H_2\end{array}/\begin{array}{c}\hfill C\mathrm{O}\end{array}=1.5$	Pellet, $\mathrm{d}=15.2 \mathrm{mm}$	R3	$117.2 \mathrm{kJ/mol}$	Isothermal	[53,72]
$\begin{array}{c}\hfill {\mathrm{Fe}}_2{\mathrm{O}}_3\end{array}\to \begin{array}{c}\hfill \mathrm{Fe}\end{array}$ (Fe2O3-20% ZrO2)	1023–1173 K	$\begin{array}{c}\hfill H_2/\mathrm{Ar}\end{array}=10$–25%	Powder, $\mathrm{m}=20 \mathrm{mg}$	R3	$59 \mathrm{kJ/mol}$	Isothermal	[59]
$\begin{array}{c}\hfill {\mathrm{Fe}}_2{\mathrm{O}}_3\end{array}\to \begin{array}{c}\hfill \mathrm{Fe}\end{array}$	973–1173 K	10% $\begin{array}{c}\hfill H_2\end{array}$ + 90% $\begin{array}{c}\hfill N_2\end{array}$	Powder, $\mathrm{m}=13 \mathrm{mg}$	A2	$23.9 \mathrm{kJ/mol}$	Isothermal	[70]
$\begin{array}{c}\hfill {\mathrm{Fe}}_3{\mathrm{O}}_4\end{array}\to \begin{array}{c}\hfill \mathrm{Fe}\end{array}$	610–877 K	100% $\begin{array}{c}\hfill H_2\end{array}$	Powder, $\mathrm{m}=100 \mathrm{mg}$	A2/A3 $<690 K$ R1 $>691 K$	39–$88 \mathrm{kJ/mol}$	Isothermal	[47]
$\begin{array}{c}\hfill {\mathrm{Fe}}_3{\mathrm{O}}_4\end{array}\to \begin{array}{c}\hfill \mathrm{Fe}\end{array}$	493–956 K	10% $\begin{array}{c}\hfill H_2\end{array}$ + 90% $\begin{array}{c}\hfill N_2\end{array}$	Powder, $\mathrm{m}=100 \mathrm{mg}$	A2/A3 $<688 K$ Rx $>691 K$	36–$103 \mathrm{kJ/mol}$	Isothermal	[47]
$\begin{array}{c}\hfill {\mathrm{Fe}}_3{\mathrm{O}}_4\end{array}\to \begin{array}{c}\hfill \mathrm{Fe}\end{array}$	473–723 K	100% $\begin{array}{c}\hfill H_2\end{array}$	Powder, $\mathrm{m}=100 \mathrm{mg}$	R2/R3	$200 \mathrm{kJ/mol}$	Isothermal	[73]
	523–663 K	100% $\begin{array}{c}\hfill H_2\end{array}$	Powder, $\mathrm{m}=100 \mathrm{mg}$	R2/R3	$71 \mathrm{kJ/mol}$	Isothermal	[73]
	>663 K	100% $\begin{array}{c}\hfill H_2\end{array}$	Powder, $\mathrm{m}=100 \mathrm{mg}$	A2	$44 \mathrm{kJ/mol}$	Isothermal	[73]
$\begin{array}{c}\hfill {\mathrm{Fe}}_2{\mathrm{O}}_3\end{array}\to \begin{array}{c}\hfill \mathrm{Fe}\end{array}$	973–1373 K	100% $\begin{array}{c}\hfill H_2\end{array}$	Cyl., $\mathrm{dxh}=9.8\times 11.1 {mm}^2$	-	$53.5 \mathrm{kJ/mol}$	Isothermal	[74]
$\begin{array}{c}\hfill {\mathrm{Fe}}_2{\mathrm{O}}_3\end{array}\to \begin{array}{c}\hfill \mathrm{Fe}\end{array}$	1123–1323 K	100% $\begin{array}{c}\hfill H_2\end{array}$	Pellet, $\mathrm{d}=10$–$15 \mathrm{mm}$	R3	$122 \mathrm{kJ/mol}$	Isothermal	[60]
$\begin{array}{c}\hfill {\mathrm{Fe}}_3{\mathrm{O}}_4\end{array}\to \begin{array}{c}\hfill \mathrm{Fe}\end{array}$	603–1003 K	5% $\begin{array}{c}\hfill H_2\end{array}$/$\begin{array}{c}\hfill N_2\end{array}$	Powder, $\mathrm{m}=10$–$15 \mathrm{mg}$	A2	$70.4 \mathrm{kJ/mol}$	Nonisothermal	[75]
$\begin{array}{c}\hfill {\mathrm{Fe}}_3{\mathrm{O}}_4\end{array}\to \begin{array}{c}\hfill \mathrm{Fe}\end{array}$	703–953 K	83% $\begin{array}{c}\hfill H_2\end{array}$/Ar	Powder, $\mathrm{m}=100 \mathrm{mg}$	-	45–$55 \mathrm{kJ/mol}$	Nonisothermal	[76]
$\begin{array}{c}\hfill {\mathrm{Fe}}_3{\mathrm{O}}_4\end{array}\to \begin{array}{c}\hfill \mathrm{Fe}\end{array}$	600–1000 K	10% $\begin{array}{c}\hfill H_2\end{array}$/Ar	Powder, $\mathrm{m}=100 \mathrm{mg}$	-	93.2–$103 \mathrm{kJ/mol}$	Nonisothermal	[77]
$\begin{array}{c}\hfill {\mathrm{Fe}}_3{\mathrm{O}}_4\end{array}\to \begin{array}{c}\hfill \mathrm{FeO}\end{array}$	450–1000 K	10% $\begin{array}{c}\hfill H2\end{array}$/Ar	Powder, $\mathrm{m}=100 \mathrm{mg}$	-	$77.3 \mathrm{kJ/mol}$	Nonisothermal	[78]
$\begin{array}{c}\hfill \mathrm{FeO}\end{array}\to \begin{array}{c}\hfill \mathrm{Fe}\end{array}$	450–1000 K	10% H2/Ar	Powder, $\mathrm{m}=100 \mathrm{mg}$	-	$85.7 \mathrm{kJ/mol}$	Nonisothermal	[78]
$\begin{array}{c}\hfill {\mathrm{Fe}}_3{\mathrm{O}}_4\end{array}\to \begin{array}{c}\hfill \mathrm{FeO}\end{array}$	973–1273 K	100% H2	Powder, $\mathrm{m}=200 \mathrm{mg}$	-	$42 \mathrm{kJ/mol}$	Isothermal	[79]
$\begin{array}{c}\hfill \mathrm{FeO}\end{array}\to \begin{array}{c}\hfill \mathrm{Fe}\end{array}$	973–1273 K	100% H2	Powder, $\mathrm{m}=200 \mathrm{mg}$	-	$55 \mathrm{kJ/mol}$	Isothermal	[79]
$\begin{array}{c}\hfill {\mathrm{Fe}}_3{\mathrm{O}}_4\end{array}\to \begin{array}{c}\hfill \mathrm{FeO}\end{array}$	1533–1573 K	$\begin{array}{c}\hfill H_2\end{array}/\begin{array}{c}\hfill C\mathrm{O}\end{array}=1$	Pellet, $\mathrm{d}=14 \mathrm{mm}$	F1	$54.1 \mathrm{kJ/mol}$	Isothermal	[80]
$\begin{array}{c}\hfill \mathrm{FeO}\end{array}\to \begin{array}{c}\hfill \mathrm{Fe}\end{array}$	1533–1573 K	$\begin{array}{c}\hfill H_2\end{array}/\begin{array}{c}\hfill C\mathrm{O}\end{array}=1$	Pellet, $\mathrm{d}=14 \mathrm{mm}$	F1	$65.4 \mathrm{kJ/mol}$	Isothermal	[80]
$\begin{array}{c}\hfill \mathrm{Fe}\end{array}\to \begin{array}{c}\hfill {\mathrm{Fe}}_3{\mathrm{O}}_4\end{array}$ (Fe2O3-20% ZrO2)	1023–1173 K	H2O/Ar = 5–20%	Powder, $\mathrm{m}=20 \mathrm{mg}$	R3	$27 \mathrm{kJ/mol}$	Isothermal	[59]
$\begin{array}{c}\hfill \mathrm{Fe}\end{array}\to \begin{array}{c}\hfill {\mathrm{Fe}}_3{\mathrm{O}}_4\end{array}$ (Fe2O3)	633–713 K	H2O/Ar = 1:4	Powder, $\mathrm{m}=1500 \mathrm{mg}$	D3	55.5–$65.3 \mathrm{kJ/mol}$	Isothermal	[81]
$\begin{array}{c}\hfill \mathrm{Fe}\end{array}\to \begin{array}{c}\hfill {\mathrm{Fe}}_3{\mathrm{O}}_4\end{array}$ (Fe2O3-8% Mo)	633–713 K	H2O/Ar = 1:4	Powder, $\mathrm{m}=1500 \mathrm{mg}$	D3	36.2–$45.2 \mathrm{kJ/mol}$	Isothermal	[81]
$\begin{array}{c}\hfill \mathrm{Fe}\end{array}\to \begin{array}{c}\hfill {\mathrm{Fe}}_3{\mathrm{O}}_4\end{array}$	1023–1173 K	$\begin{array}{c}\hfill H_2\mathrm{O}\end{array}=10 {\mathrm{cm}}^3\left(\mathrm{STP}\right)/\min$,	Powder, $\mathrm{m}=1.0 \mathrm{mg}$	D3	$11 \mathrm{kJ/mol}$	Isothermal	[57]
		$\begin{array}{c}\hfill N_2\end{array}=20 {cm}^3\left(\mathrm{STP}\right)/\min$
$\begin{array}{c}\hfill \mathrm{Fe}\end{array}\to \begin{array}{c}\hfill {\mathrm{Fe}}_3{\mathrm{O}}_4\end{array}$ (Fe2O3)	773–1173 K	$\begin{array}{c}\hfill H_2\mathrm{O}\end{array}$	Powder, $\mathrm{m}=7$–8 mg	D3	$77.9 \mathrm{kJ/mol}$	Isothermal	[82]

. In the following, the individual results of the literature review are discussed, using the symbols from Table 1.

The literature review reveals that although some researchers have chosen nucleation models (A2) as the control mechanism, chemical reaction at the interfaces (R2, R3) and reaction order-based models (F1) are also commonly applied models within the related studies. Furthermore, some researchers have shown that the rate-controlling step can be the reduction of wuestite (FeO) to iron [53,70]. The formation of a dense iron shell on the surface of the wuestite layer makes the diffusion of hydrogen through this shell more difficult. Above $576 °\mathrm{C}$, the reduction of magnetit Fe₃O₄ by hydrogen H₂ proceeds in two steps, according to the Baur–Glaessner diagramm (Figure 6) and the following simplified reduction equations [71]:

$${\mathrm{Fe}}_3{\mathrm{O}}_4+H_2\to 3FeO+H_2\mathrm{O}$$

(43)

$$\mathrm{FeO}+H_2\to \mathrm{Fe}+H_2\mathrm{O}$$

(44)

Wuestite is an intermediate product in the reduction of magnetit. Below $576 °\mathrm{C}$, wuestite is unstable when magnetite is reduced by hydrogen, as shown in Figure 6. In this range, the reduction takes place according to Equation (2).

One of the most influencing parameters on the reaction kinetics is the activation energy $E_{\mathrm{a}}$. It is essential for determining, for example, the reactor dimension and the energy consumption to sustain the reaction, e.g., at a constant temperature. Based on the literature review for the reduction of iron oxide (see Table 5), $E_{\mathrm{a}}$ as well as the dominant kinetic models $f\left(\alpha \right)$ depends on the starting raw material, the type and concentration of the reaction gas, the reaction temperature, impurities as well as the physical form of the reacting solid material. Because of these numerous influencing parameters, activation energies in the range from 36 to $200 \mathrm{kJ/mol}$ have been reported in the literature.

Even small changes in the sample structure or in the process parameters have a significant effect on the kinetics data. For this reason, the next sections describe and discuss the kinetic mechanisms as well as the activation energy for the investigated storage material upon changing the reaction temperature and the reacting gas concentration.

4.1.1. Kinetic Mechanism and Effect of the Process Temperature

As shown in Figure 5, the first step towards the estimation of the kinetic parameters is the selection of the kinetic mechansim from Table 1. To this aim, the experimental data have been normalised and plotted as shown in Figure 7, showing the isothermal reduction curves for $600 °\mathrm{C}$, $700 °\mathrm{C}$ and $800 °\mathrm{C}$ over time (Figure 7a) and the time derivative of the reduction fractions over the reduction fraction (Figure 7b) of the sample mass. The presented results in Figure 7 belong exemplarily to the reactant gas concentration ${\mathrm{ro}}_3$ according to Table 4.

Increasing the reaction temperature of the sample from 600 to $700 °\mathrm{C}$ results in a decrease in the time needed for the reduction fraction to reach ${\alpha }_{\mathrm{Re}}=0.8$ from 134 to $43 \min$ (68%). Increasing the temperature further to $800 °\mathrm{C}$ results in an additional acceleration of the reduction reaction by 47% (from 43 to $23 \min$). It can be concluded that the reduction reaction becomes faster with increasing the reaction temperature, which is common for such an endothermic reaction [45,83,84]. Furthermore, the curves in Figure 7b show a decelerating behavior, which indicates a phase boundary controlled reaction model. The high gradient, at the beginning, can be attributed to a short-term dominant surface diffusion, as illustrated in Figure 3 (Deceleratory).

First assumption from the graphical interpretation of Figure 7b can be confirmed by comparing the relevant reaction mechanism model from Table 1 in Figure 8. With the geometrical contracting sphere model R3 (see Figure 8c), the discrepancy between measured to fitted curves at a reactant flow rate of $250 \mathrm{mL/min}$ with gas concentration ro${}_3$ (see Table 4) and reaction temperature of $700 °\mathrm{C}$ amounts to AD${}_{\mathrm{Re},700 °\mathrm{C},\mathrm{R}3}$≤ 1.9%. Furthermore, the AD upon applying the A1.5 model decreases to less than AD${}_{\mathrm{Re},700 °\mathrm{C},\mathrm{A}1.5}<$ 1% from $56 \min$ and a conversion degree of $\alpha >0.919$. This result might indicate that nucleation growth can be the dominating process from this point on. As the conversion rate is already quite high, such an effect might not be important for a later practical application. However, from the scientific point of view, more dedicated investigations are needed to find out whether a hybridisation between different models might be leading to a much higher modelling accuracy. In addition, the model R2 reveals also a moderate deviation from AD${}_{\mathrm{Re},700 °\mathrm{C},\mathrm{R}2}\le$ 2.4% over the reaction time. Applying the other reaction models F1, A1.5 (Time $\le 56 \min$) and D3 results in remarkably increased deviations between experimental and theoretical curves as depicted in Figure 8.

In addition, an examination of the time derivative of the reduction fraction over the reduction fraction in Figure 9b also shows, above a conversion rate of $\alpha =0.05$, a maximum deviation of AD${}_{\mathrm{Re},700 °\mathrm{C},\mathrm{R}{3}^{\prime }}\le 5\%$, which proves that R3 or the geometrical contracting sphere model is perfectly suitable to describe the experimentally obtained hydrogen reduction kinetics at $700 °\mathrm{C}$.

Figure 9 and Figure 10 confirm the excellent suitability of the reaction mechanism model R3 also for reaction temperatures $600 °\mathrm{C}$ and $800 °\mathrm{C}$. With the model R3, the difference between the experimental curves and the fitted ones amounts to AD${}_{\mathrm{Re},600 °\mathrm{C},\mathrm{R}3}<$ 2% and AD${}_{\mathrm{Re},800 °\mathrm{C},\mathrm{R}3}<$ 2.5%.

4.1.2. Effect of the Process Gas Composition

Figure 11a depicts the relative reduction fractions of the investigated pellet with different H₂/N₂ concentrations (cf. Table 4) at $700 °\mathrm{C}$. The time needed to achieve a reduction of 80% of the maximum possible mass change (${\alpha }_{\mathrm{Re}}=0.8$) decreases to less than 72%, from $160 \min$ down to $45 \min$ as the H₂ concentration in the reactive gas mixture increases from ro${}_1$ to ro${}_3$. This is also confirmed by the study of [84]. Furthermore, the time derivative of the reduction fraction over the reduction fraction in Figure 11b shows the same decelerating behavior as Figure 11a, confirming the most suitability of a geometric contraction model.

Applying the geometrical contracting model R3 (see Figure 12a), the deviation between measured and fitted curves with gas concentration ro${}_2$ and reaction temperature of $700 °\mathrm{C}$ amounts to AD${}_{\mathrm{Re},{\mathrm{ro}}_2,\mathrm{R}3}\le$ 1.4%.

The investigation of temporal change in the reduction fraction over the reduction fraction in Figure 12b also shows a very good agreement between measured and calculated data. From a conversion rate of $\alpha =0.05$, the discrepancy is AD${}_{\mathrm{Re},{\mathrm{ro}}_2,\mathrm{R}{3}^{\prime }}\le$ 5%.

Figure 13 confirms the very good fit of the reaction mechanism model R3 to describe the experimental reduction kinetic data for gas concentration ro${}_1$–ro${}_3$ together with Figure 12. With the model R3, the difference between the experimental curves and the fitted ones amounts to AD${}_{\mathrm{Re},{\mathrm{ro}}_1,R3}<$ 1.3% and AD${}_{\mathrm{Re},{\mathrm{ro}}_3,\mathrm{R}3}<$ 1.9%. The optimised rate constants $k_{\mathrm{A},n}$ of Equation (40) for the R3 model are given on the top line of Figure 12a and Figure 13. It can be concluded that the higher the gas concentration, the higher the reduction rate constant $k_{\mathrm{A},\mathrm{Re}}$.

4.2. Oxidation

In contrast to the reduction reaction, there exist considerably less studies in the literature on the steam oxidation of iron for hydrogen production (see Table 5). This is mainly due to the fact that a large number of iron oxide reduction studies concentrate on the reduction of iron oxide in iron production processes. However, in the last few years, the number of publications dealing with hydrogen production by means of iron oxidation has increased [57,59,81,82]. Most of these research articles assumed that the oxidation is controlled by diffusion through the oxide layer, which is described by the diffusion model D3 (see Table 1). Kang [59], on the other hand, has shown that the control mechanism of the oxidation reaction can be determined by the chemical reaction at the interface (R3) upon studying a powdered calcined Fe₂O₃-20%ZrO₂ sample.

In general, the oxidation of iron is the reversible process of Equations (2), (43) and (44). If the reaction temperature is above $576 °\mathrm{C}$ (the eutectoid point of the Fe-O system, see Figure 6), multilayer growth from pure iron to wuestite to magnetite takes place. Below $576 °\mathrm{C}$, the wuestite phase is unstable, and the oxidation of iron produces directly magnetite, as already described for the reduction reaction in Section 4.1. For a detailed description of the individual oxidation steps and the associated valence ion migration, references [85,86] are highly recommended.

In addition, for the oxidation reaction, a relatively wide range for the oxidation activation energy (from 11 to $77.9 \mathrm{kJ/mol}$) exist in the literature as presented in the lower part of Table 5. In order to have a precise simulation model for subsequent research tasks with the investigated storage material, a set of isothermal oxidation experiments was conducted at temperatures $600 °\mathrm{C}$, $700 °\mathrm{C}$ and $800 °\mathrm{C}$. In addition, the reacting gas concentration have been varied from ro${}_1$ to ro${}_3$ during the test campaign (see Table 4).

4.2.1. Kinetic Mechanism and Effect of the Process Temperature

As before, the kinetic mechanism formalism must be defined first in order to subsequently determine the kinetic parameters. For this purpose, the measured oxidation reaction data were normalised and plotted as shown in Figure 14. Here, the experimentally measured isothermal oxidation curves for $600 °\mathrm{C}$, $700 °\mathrm{C}$ and $800 °\mathrm{C}$ over time are plotted (Figure 14a) and the time derivative of the oxidation fractions over time versus oxidation fraction (Figure 14b) for the reactant gas concentration ${\mathrm{ro}}_3$ are shown in Table 4.

The curves in Figure 14a demonstrate an acceleration of the kinetic rate with the increase of the reaction temperature, which is a similiar result to those reported in [57,82]. Increasing the temperature of the sample from 600 to $700 °\mathrm{C}$ leads to a reduction in the time needed for the oxidation fraction to reach ${\alpha }_{\mathrm{Ox}}=0.8$ from 41 to $35 \min$ (15%). Increasing the temperature further to $800 °\mathrm{C}$ leads to an additional acceleration of the oxidation reaction by 16% (time needed to reach a normalised conversion rate of 80% decreases from 32 to $27 \min$).

Furthermore, the derivation of the measured curves represented in Figure 14b show an initial sigmoidal behavior between ${\alpha }_{\mathrm{Ox}}=0$ and ${\alpha }_{\mathrm{Ox}}=0.18$. After attaining a maximum value, a deceleratory shape is observable, which indicates a phase boundary controlled reaction model (see Figure 3).

Different reaction mechanism models out of Table 1, namely, F1, R2, R3, D3 and A1.5, have been tested on their suitability to correlate the measured oxidation kinetic data, and the results are depicted in Figure 15. The presented curves belong to the isothermal oxidation temperature of $700 °\mathrm{C}$ and the reactant flow rate of $250 \mathrm{mL/min}$ at the gas concentration ro${}_3$ (see Table 4).

Sigmoidal behavior up to ${\alpha }_{\mathrm{Ox}}=0.18$ indicates a nucleation model, which is also consistent with the small deviation from AD${}_{\mathrm{Ox},700 °\mathrm{C},\mathrm{A}1.5}\le$ 1% of the Avrami model A1.5 in Figure 15e. The deceleratory slope from ${\alpha }_{\mathrm{Ox}}=0.18$ to ${\alpha }_{\mathrm{Ox}}=1$ is again most accurately represented by the geometrical contracting model R3 in Figure 15c. Here, the deviation between measured and fitted amounts to AD${}_{\mathrm{Ox},700 °\mathrm{C},\mathrm{R}3}\le$ 1.8%. Other models, such as reaction-based model F1 (Figure 15a) showed a deviation in the range AD${}_{\mathrm{Ox},700 °\mathrm{C},\mathrm{F}1}<$ 11%. The geometrical contracting model R2 (Figure 15b) is associated with a deviation in the range of AD${}_{\mathrm{Ox},700 °\mathrm{C},\mathrm{R}2}<$ 6%. The worst correlation results have been obtained upon applying the diffusion model D3 (Figure 15d) showing a deviation of AD${}_{\mathrm{Ox},700 °\mathrm{C},\mathrm{D}3}\le$ 30%. In Figure 15f, the optimised rate constants $k_{\mathrm{A},n}$ of Equation (40) for all the tested reaction mechanism models are presented.

Figure 16 demonstrates the very good fitting of the reaction mechanism model R3 to correlate the experimental oxidation kinetic data at the temperatures $600 °\mathrm{C}$, $700 °\mathrm{C}$ and $800 °\mathrm{C}$. The correlation has been carried out above the conversion rate of ${\alpha }_{\mathrm{Ox}}>0.01$, which adjusts $k_{\mathrm{A},\mathrm{Ox},n}$. A detailed explanation for the reason behind this truncation is given in the next subsection. The deviation between the experimental curves and the fitted ones is AD${}_{\mathrm{Ox},600 °\mathrm{C},\mathrm{R}3}<$ 3.2%, AD${}_{\mathrm{Ox},700 °\mathrm{C},\mathrm{R}3}<$ 2.0% and AD${}_{\mathrm{Ox},800 °\mathrm{C},\mathrm{R}3}\le$ 1.2%.

4.2.2. Effect of the Process Gas Composition

Figure 17a illustrates the oxidation conversion rates of the investigated pellet with different H₂O/N₂ concentrations (cf. Table 4) at $700 °\mathrm{C}$. The time needed to achieve a conversion rate of ${\alpha }_{\mathrm{Ox}}=0.8$ decreases to less than 66% from $82 \min$ down to $28 \min$ as the H₂O concentration in the reactive gas mixture increases from ro${}_1$ to ro${}_3$.

Furthermore, the time derivatives of the oxidation fractions over the oxidation fraction in Figure 17b exhibit a non-consistent course in the range of ${\alpha }_{\mathrm{Ox}}\le 0.13$. For this purpose, the additional extract in Figure 17a is presented. It is obvious that a changing course of the temporal conversion rate variation exists for the operating condition ro${}_3$, if compared with the other two operating conditions (ro${}_1$ and ro${}_2$), over the first $5.5 \min$. Since the cubic starting slope of the conversion rate changes to a linear behavior with falling concentration, it is assumed that the non-consistency of ro${}_3$, at the begining, appears because of the dead time between the switchover of the MFC from pure N₂ flow for the inertisation of the measuring cell to the flow time point of the H₂O/N₂ mixture. The steam generator needs, obviously, a couple of minutes until a continuous steam mass flow can be realised. The subsequent slope (${\alpha }_{\mathrm{Ox}}>0.18$) indicates a decelerating behavior as shown in Figure 14b, confirming the geometric contraction model. For the fitting purpose, the truncation of the data before ($\alpha =0.01$) for this operating condition was enough to obtain a very good consistency with the other operating conditions.

Using the geometric contraction model R3 (see Figure 18), the maximum deviation between measured and fitted values at a reaction temperature of $700 °\mathrm{C}$ is AD${}_{\mathrm{Ox},\mathrm{j}}\le$ 2.7% (AD${}_{\mathrm{Ox},{\mathrm{ro}}_1,\mathrm{R}3}\le$ 2.7%, AD${}_{\mathrm{Ox},{\mathrm{ro}}_2,\mathrm{R}3}\le$ 1.8%, AD${}_{\mathrm{Ox},{\mathrm{ro}}_3,\mathrm{R}3}\le$ 2%). Therefore, the functionality of the reaction mechanism model R3 for gas concentration ro${}_1$–ro${}_3$ is confirmed. Furthermore, the higher the gas concentration, the higher the oxidation rate $k_{\mathrm{A}},\mathrm{Ox},\mathrm{n}$, as depicted on the header line of each separate figure.

4.3. Activation Energy and Kinetic Velocity Rate

For developing a precise thermal management system of a redox reactor in the following development steps, the activation energy shall be described as accurately as possible as a function of its influencing parameters. Among other parameters, the activation energy $E_{\mathrm{a}}$ depends on the chemical composition of the solid material, its physical properties and process temperature. All these parameters have an influence on the rate-limiting step, and therefore, different rate-limiting steps lead to different activation energies. The estimated values of the activation energy related to the hydrogen reduction of iron oxides and steam oxidation of iron in various studies are listed in Table 5. Indeed, the listed values for the activation energy exhibit a wide range for the reduction reaction, from 36 to $200 \mathrm{kJ/mol}$, as well as for the oxidation reaction, 11 to $77.9 \mathrm{kJ/mol}$. In addition, $E_{\mathrm{a}}$ for both reactions, reduction and oxidation, decreases with increasing the temperature. The oxidation process shows also a higher reactivity with increasing the steam concentration.

It can, therefore, be concluded that the activation energy cannot be selected from the literature unless the material and the operating conditions are exactly similar. In other words, for the introduced iron oxide pellet in this work, the activation energy shall be described as a function of its influencing parameters based on the obtained experimental kinetic data. One of the most common methodologies to describe the activation energy is its formulation as an Arrhenius equation. Arrhenius published the followed equation in the year 1889 based on empirical considerations for the calculation of the activation energy [87]:

$$k_{\mathrm{A}}=A{e}^{-\frac{E_{\mathrm{a}}}{RT}}$$

(45)

The temperature dependency is given by the exponential term and the pre-exponential factor A is a constant. It is possible to determine the activation energy from a series of kinetic measurements by plotting the logarithm of the rate constant against the reciprocal of the absolute temperature, as shown by the transformation of the Equation (35) to

$$ln\left(k_{\mathrm{A}}\right)=ln\left(A\right)-\frac{E_{\mathrm{a}}}{R}\frac{1}{T}$$

(46)

The slope of such a straight line is proportional to the activation energy and the intercept to the pre-factor (A). The following equation can be simply derived from the past equation to estimate the activation energy.

$$E_{\mathrm{a}}=-R\frac{\partial ln\left(k_{\mathrm{A}}\right)}{\partial (1/T)}=R{T}^2\frac{\partial ln\left(k_{\mathrm{A}}\right)}{\partial T}$$

(47)

Combined with Equation (41), the activation energy graphs for the reduction reaction (see Figure 19) and the oxidation reaction (see Figure 20) were calculated. Figure 19a and Figure 20a illustrate the linearised descriptions according to Equation (46) for the reduction and oxidation reactions, respectively, while Figure 19b and Figure 20b depict the general Arrhenius form according to Equation (45). From the slopes of the linearised graphs, the activation energy for reduction is estimated to $E_{\mathrm{a},\mathrm{Re}}$ = 56,858 J/mol and for oxidation to $E_{\mathrm{a},\mathrm{Ox}}$ = 16,001 J/mol. The estimated $E_{\mathrm{a}}$ values can be utilised to describe the temperature-dependency of the rate constant for both reduction and oxidation with a coefficient of determination of $R^2=0.98336$ and $R^2=0.9858$, respectively.

Compared to reduction, the temperature has a 71% smaller effect on the reaction kinetics of oxidation. This is attributed to the lower activation energy of oxidation. Furthermore, both reaction rates increase with increasing the temperature. In terms of the reduction reaction, the obtained temperature dependency is consistent with the endothermic reaction mechanism. On the other hand, the temperature dependency of the oxidation reaction is traced back to the fact that enhanced diffusion processes have a larger impact on the reaction rate than the inhibited exothermic reaction mechanism due to the higher temperatures.

In addition to the discussed temperature dependency of the velocity rate $k_{\mathrm{A}}$, the dependency on the reaction gas concentration C can also be described in the general Arrhenius equation. To this aim, the concentration C raised to an empirical exponent parameter n is added as a multiplier to the typical Arrhenius equation to describe both temperature and concentration of the respective reactive gas component in the gas phase on the reaction rate [58,59,60]:

$$k_{\mathrm{A}}(C,T)=A{e}^{-\frac{E_{\mathrm{a}}}{RT}}{C}^n$$

(48)

In combination with the calculated activation energies and the non-linear least square fitting methodology (see Figure 5), the parameter n was determined to be $n=0.75$ for both reduction and oxidation. With a coefficient of determination of $R^2=0.9211$ and $R^2=0.958$, Figure 21a,b illustrate a good accuracy with the determined parameter n for the concentration-dependent description of the velocity rate at, for instance, a reaction temperature of $700 °\mathrm{C}$. Furthermore, according to Figure 21a,b, both reaction rates increase with increasing the concentration of the respective reactive gas. In case of oxidation, the incremental rate of change with rising concentration is almost twice as high as in the case of reduction. Indeed, more experimental runs are needed at other reacting gas concentrations and temperatures in order to further enhance the fitting accuracy ($R^2$) of the model. This will be the subject of our near future activities and publications.

In summary, based on the parameters for $E_{\mathrm{a}}$, A, and n determined from the redox experiments described here, a temperature- and concentration-dependent description for the velocity rate constant can be expressed as follows:

$$\begin{array}{c}{k}_{\mathrm{A},\mathrm{Re}}(C,T)=0.1813 s^{-1}·C_{{\mathrm{H}}_2\mathrm{O},\mathrm{in}}^{0.75}·e^{-\frac{56.9 \mathrm{kJ/mol}}{RT}}\end{array}$$

(49)

$$\begin{array}{c}{k}_{\mathrm{A},\mathrm{Ox}}(C,T)=0.0023 s^{-1}·C_{{\mathrm{H}}_2\mathrm{O},\mathrm{in}}^{0.75}·e^{-\frac{16.0\mathrm{kJ/mol}}{RT}}\end{array}$$

(50)

Figure 22a,b present the calculated velocity rate constants of Equations (49) and (50) as a function of temperature T and gas concentration C. Compared to oxidation, the temperature has a smaller effect on the velocity rate constant of reduction, which is due to the higher activation energy. This explains also the lower temperature gradients in Figure 22a compared to Figure 22b.

Based on the identified best-suited kinetic model (R3) and the estimated parameters $E_{\mathrm{a}}$, A, and n, the following expressions are suggested to describe the temperature- and concentration-dependent conversion rate of both reduction and oxidation reactions on the investigated iron oxide pellet:

$$\begin{array}{c}\frac{d{\alpha }_{\mathrm{Re}}(C,T)}{dt}=0.5439 s^{-1}·C_{{\mathrm{H}}_2\mathrm{O},\mathrm{in}}^{0.75}·e^{-\frac{56.9 \mathrm{kJ/mol}}{RT}}·{(1-\alpha )}^{2/3}\end{array}$$

(51)

$$\begin{array}{c}\frac{d{\alpha }_{\mathrm{Ox}}(C,T)}{dt}=0.0069 s^{-1}·C_{{\mathrm{H}}_2\mathrm{O},\mathrm{in}}^{0.75}·e^{-\frac{16.0 \mathrm{kJ/mol}}{RT}}·{(1-\alpha )}^{2/3}\end{array}$$

(52)

5. Conclusions

In this work, the kinetics of both reduction and oxidation (redox) reactions on an iron oxide/iron pellet have been experimentally investigated under different isothermal conditions and different concentrations of the respective reacting gas. In addition, a comprehensive literature review has been conducted on the existing kinetic models and their suitability to describe the obtained results. It can be concluded that the so-called shrinking core model is best suited to describe the measured kinetic data. Additionally, the literature review showed a wide range for the activation energies of both redox reactions. Furthermore, the reported kinetic parameters are quite sensitive to small changes in the sample structure as well as in the process parameters. For this reason, the kinetic mechanisms and the activation energy for the investigated storage material were determined under different reaction temperatures and process gas concentrations.

Furthermore, it turned out that by applying the geometrical contracting sphere reaction model (R3), both measured concentration- and temperature-dependent reduction and oxidation rates could be reproduced analytically with a maximum deviation of less than 5%. Higher concentrations of the relevant reactive gas component within the input gas stream as well as elevated temperatures resulted in an accelerated conversion and, accordingly, a shorter time period to accomplish the reaction. Compared to the oxidation reaction, the temperature showed a greater effect on the reduction kinetics, which is attributed to the higher reduction activation energy ($E_{\mathrm{a},\mathrm{Re}}=56.9 \mathrm{kJ/mol}$ versus $E_{\mathrm{a},\mathrm{Ox}}=16.0 \mathrm{kJ/mol}$). Oxidation showed, however, a nearly twofold increase in the reaction rate constant with increasing concentration if compared to the reduction reaction.

With the estimated activation energies, the temperature dependency of the rate constants for both reduction and oxidation could be very well described with an Arrhenius equation with a coefficient of determination of $R^2=0.98336$ and $R^2=0.9858$, respectively. To describe the concentration C dependence of the rate constants, an empirical exponent has been introduced and estimated to $n=0.75$, with which the velocity rate $k_{\mathrm{A}}$ can be well described with a coefficient of determination for the reduction of $R^2=0.913$ and for the oxidation of $R^2=0.958$.

The obtained results pave the way to implement the obtained correlations to model both isothermal and non-isothermal behaviours of single- and multi-pellet redox reactors for designing and optimising the hydrogen storage process. These development results will be the subject of our next research.