A. S. Portnyagin*ab,
A. P. Golikova,
V. A. Drozdc and
V. A. Avramenkoab
aDepartment of Sorption Processes, Institute of Chemistry, Far Eastern Branch of Russian Academy of Sciences, pr. 100-letiya Vladivostoka, Vladivostok, Russia. E-mail: arsuha@gmail.com
bSchool of Natural Sciences, Far Eastern Federal University, Sukhanova str 8., Vladivostok, Russia
cScientific Educational Center of Nanotechnology, Far Eastern Federal University, Sukhanova str 8., Vladivostok, Russia
First published on 16th January 2018
To date, kinetic computations have been carried out efficiently for a great variety of physico-chemical processes including crystallization, melting and solid–solid transitions. However, appropriate methods for the kinetic analysis of chemical reactions, especially multi-staged reactions, are currently lacking. Here we report on an alternative way of treating temperature-programmed reaction data using the reduction of iron(III) oxide as an example. The main principle in the suggested approach is to take into account every stage of the studied process, resulting in a system of kinetic differential equations. Kinetic parameters (activation energy and preexponential factors) are optimized for each of the stages, and cubic splines are used to approximate the conversion functions that reflect changes in reaction-specific surface area throughout the process. The applicability of the suggested method has been tested on temperature-programmed reduction (TPR) data for iron(III) oxide samples produced from the original Fe2O3 powder by annealing it at 600, 700 and 800 °C. Results of kinetic analysis obtained at different temperature regimes demonstrate the good stability and performance of the method. Peculiarities of iron(III) oxide reduction have been revealed, depending on the stage and heating rate. The influence of material morphology on the reduction kinetics has been assessed by comparing preexponential factors corresponding to the first reduction stage. This approach allows a comparison of the structural characteristics of the materials based on the kinetic analysis of the TPR data. Using optimized conversion functions, the initial particle size distribution has been reproduced. Theoretically found particle size distribution was found to correlate well with the experimental distribution obtained via laser diffraction.
Despite the wide scope of TPR applications in investigational practice, analysis of experimental results is carried out far more superficially. As usual, conclusions are based on such facile assessments as analysis of peak shape alterations or changes in the location of the peak maxima of the TPR curve.13–15 However, more detailed information about the underlying processes of temperature-programmed reactions can be obtained by a quantitative description of the TPR processes conducted via kinetic analysis.16
There are two main groups of method used to process non-isothermal experimental data: isoconversional and model fitting methods. These methods are commonly based on eqn (1):
(1) |
The basic assumption of isoconversional (model-free) methods is that the reaction rate depends only on temperature at the same degree of conversion. In order to implement this assumption, it is necessary to have experimental data recorded at several heating rates. Processing non-isothermal experimental data is carried out by means of various model-free methods.17 In the method suggested by Friedman,18 the logarithmic representation of rate eqn (1) is used:
(2) |
Then, for a chosen α, the plot ln(βdα/dT) vs. T−1 obtained from the TPR curves recorded at different heating rates is built up. If it is a straight line, its slope allows evaluation of the activation energy.
The Flynn–Wall–Osawa method,19 based on the integral form of eqn (1), includes Doyle's approximation for temperature integral:20
(3) |
The procedure for evaluating activation energies is similar to the procedure for the Friedman method.
If the process is characterized by one rate-limiting step or several separated in time, then it is possible to obtain a unique set of kinetic parameters. However, in the case of a blend of oxides, in which reduction is not a separate process, or oxide reacting through multiple intermediate stages (e.g. reduction of iron oxide), kinetic analysis performed by isoconversional methods gives us variable (apparent) Ea, which depends on α. Although it is useful to use isoconversional kinetics to reveal the process mechanism, evaluation of the apparent activation energy and conversion function is complicated due to the absence of sound theory behind chemical solid–gas reactions. Even though there are a lot of mathematical expressions for f(α), each of which refers to a certain type of topochemical mechanism, it is difficult to choose between them due to the fact there can be several mechanisms acting when the wide temperature region is studied. Also, several kinetic models can fit well and be equally accurate but correspond to different activation energies and preexponential factors. Further, in the case of the multistage process mentioned above, each stage can be characterized by an individual kinetic triplet (Ea, A, f(α)), which cannot be addressed by eqn (1).
In the case of a multistage process, it is advisable to implement methods of model fitting which imply a numerical solution to the system of differential equations, each of which corresponds to one of the stages.21 Thus, the TPR curve is simulated while the kinetic parameters are adjusted to fit the experimental curve. However, because of the correlation between kinetic parameters, results of model fitting strongly depend on the initial approximation.22,23 Consequently, if no additional assumptions about the process are made, it is impossible to achieve the unique set of kinetic parameters via model fitting methods.
Here we present an attempt to overcome the aforementioned disadvantages of existing methods of kinetic analysis. The aim of this study was to create a universal approach for the kinetic analysis of non-isothermal solid–gas reaction experimental data, based on minimal model assumptions. The key point of the suggested method is to reproduce functions of “active” specific surface area S(α) (where α is conversion degree) from the set of TPR experimental curves. The suggested method allows separating temperature (Arrhenius) dependence and temperature independent (specific surface area changes) parts of the TPR curve and, therefore, obtaining the possibility of determining stable (mathematically) activation energies of separate stages in a multistage TPR process. TPR experimental data of iron(III) oxide were selected to test the suggested method.
The Brunauer–Emmett–Teller (BET) surface area of the iron oxides was determined on a Autosorb IQ (Quantachrome, USA) device. Scanning electron microscope (SEM) images were taken using a Carl Zeiss Crossbeam 1540xb (Germany) instrument. Particle size distributions were obtained using Fritsch particle sizer Analysette 22 (Germany). TPR measurements were carried out on an automated chemisorption analyzer ChemBET Pulsar TPR/TPD (Quantachrome Instr., USA). Sample powders (35–40 mg) were loaded into the quartz sample cell. To separate water forming during oxide reduction, a liquid nitrogen cold trap was used. Prior to the experiment, samples were annealed at 350 °C in a nitrogen flow for 30 minutes to degas and remove moisture. TPR curves were recorded under different temperature programs, at a flow rate of 50 mL min−1, under flow of 6% H2 + N2 and ambient pressure. Purity of gases was 99.995 vol%.
In TPR experiments, small amounts of solid material and high gas flow rates are usually used. In connection to this, the following assumptions are made: (1) because reaction by-products (in the case of TPR–water) are removed from the reaction zone with gas flow, reverse reaction can be neglected; (2) curve broadening caused by longitudinal diffusion is not to be taken into account; (3) the elementary stage of the reduction reaction is first order with respect to hydrogen.
Let's consider the kinetic model of the reduction of bivalent metal with hydrogen that corresponds to the scheme MeO + H2 = Me + H2O. Hydrogen consumption over time, dτ, is defined by eqn (4):
(4) |
n = n0(1 − α) | (5) |
Substituting (5) into (4) and taking into account that α is a single-valued function of τ we get the following equation:
(6) |
Comparing eqn (6) and (1), and assuming k(T) = A*exp(−Ea/RT), we find that:
f(α) = s(α)(1 − α)Pg | (7) |
Thus, the conversion function in that case is connected with the function of specific surface area change; i.e. if f(α) is known, then the s(α) can be evaluated and vice versa. At the present time, there are plenty of mathematical forms to express f(α).24 The choice of a certain type of conversion function is based either on a priori knowledge about the geometry of the interfacial reacting surface or just on the lowest residual dispersion between the theoretical and experimental curves. However, in both cases the chosen conversion function far from always reflects the real surface geometry of the investigated system: a priori considerations can be wrong; function choice based on minimal residual dispersion also can be wrong if none of the test functions describe the surface geometry. In connection with this, we suggest another approach to the choice of optimal conversion function. The idea is simple: instead of determination of what known function exactly corresponds to the investigated system, we evaluate the optimal function of the surface area change, corresponding to the investigated TPR curves. Further, the realization of the suggested approach will be given taking into account experimental parameters (gas flow rate, hydrogen concentration in the gaseous mixture, heating rate etc.).
(8) |
(9) |
Within the small time interval that corresponds to the time needed to fill the tube with the gaseous mixture, when the temperature change and convective transfer can be neglected, the process is described by the following system of differential equations:
(10) |
Let's introduce the modeling step, which corresponds to the time interval needed to fill the reactor's volume with the gaseous mixture. For a tubular reactor we have:
(11) |
We propose that within the time interval Δτ, the temperature does not change and equals:
(12) |
After Δτ has passed, the temperature changes abruptly up to the value Ti+1 = T(τi+1) = T(τi + Δτ). On every time interval, the system of differential equations is solved numerically, and the change in hydrogen partial pressure and the composition of the reactive mixture are estimated. During the next-time step, the reaction mixture is filled with a new portion of gaseous mixture and the process of modeling is repeated. Because Δτ is small, such a scheme allows the avoidance of a numerical solution of a differential equation of convective transfer and, at the same time, provides sufficient accuracy of the modeling (taken from preliminary numerical experiments). As a result, the model TPR is built:
(13) |
Because the time step is small, the obtained model curve can be considered as continuous. To solve the inverse problem; i.e. to identify the set of parameters A*,E,λ0,…,λn−1 that allow such model curves to fit the experimental curves (with various heating rates) as best as possible, one has to minimize the following function:
(14) |
(15) |
Let's introduce new function φ+(α) as the following:
(16) |
The following function is used as a stabilizer in the present work:
(17) |
Finally, the object function (14) is transformed into:
(18) |
After we obtained an optimal set of kinetic parameters and splines of relative specific surface area we tried to evaluate the initial particle size distribution, ω, based on the assumption that initial particle size distribution defines explicitly how relative specific surface area changes throughout the process of TPR. We take into account only the initial iron oxide, Fe2O3, because its reduction is characterized by the process of contraction of the particles during TPR, while the other oxides proceed via the growth and contraction stages that cannot be evaluated and separated precisely. In order to obtain an initial particle size distribution using the optimal s(α) spline and n(α) dependence, we conducted an iterative procedure until the particle size distribution on the nth step reproduces the optimal s(α) spline. The following iterative formula was used:
(19) |
In present work, we conducted research on the applicability of the suggested method to describe the kinetics of a heterogeneous chemical reaction on the example of TPR of iron(III) oxide. We provide the results of calculations and compare them with the corresponding experimental data.
TPR of iron(III) oxide is very complicated in terms of kinetic analysis because of the multistage character of the reduction (Fe2O3 → Fe3O4 → FeO → Fe).37 Every stage of this reduction process is characterized by a separate set of kinetic parameters that makes model-free methods, based on eqn (1), inapplicable to the kinetic analysis of such a process. This is also proved by the wide range of activation energies (18–246 kJ mol−1)32 caused not only by the various compositions and morphologies of the investigated sample, but also by inaccurate kinetic analysis. In our work we used the Friedman method to obtain dependence Ea from α in order to compare it with the corresponding values evaluated via the suggested method. The obtained curves (Fig. 1) clearly demonstrate that the reduction passes through three consecutive steps as there are three regions on the Friedman plots, where Ea values remain more or less constant. The disadvantage of the Friedman method in this case is that Ea cannot be evaluated precisely for each stage due its variation with respect to the studied material or the degree of conversion. In particular, the activation energy of the first reduction step Fe2O3–Fe3O4, proceeding while the conversion degree has not reached 0.11, lies in the range from 7 to 115 kJ mol−1, therefore making evaluation of this parameter inaccurate. Such a wide confidence interval occurs because the reduction of each oxide contributes to the dα/dt term. Because reduction reactions proceed with different rates, the contribution of each reduction reaction to the dα/dt term varies when the heating rate changes, thus widening the interval of the activation energy values. This complication arises from eqn (1), which can be utilized properly only in the case of a single-stage processes. Also, it can be seen that the activation energies of the other reduction stages oscillate from one curve to another, thus suggesting that the single-rate eqn (1) cannot be applied to a multistage process.
In our study we used samples of different morphologies to assess the variation in the kinetic analysis results with respect to the studied material. To ensure that all samples have identical chemical composition, we annealed iron(III) oxide at different temperatures to obtain series of samples with different morphologies.
SEM images of the initial iron oxide powder are shown in Fig. 2. The original powder consists of spherical particles aggregated from flakes of Fe2O3. This sample possesses a highly-developed surface structure, proved by the highest BET surface area among all studied materials (25 m2 g−1; Table 1).
Sample | EFe2O3/kJ mol−1 | EFe3O4/kJ mol−1 | EFeO/kJ mol−1 | AFe2O3/s−1 | Sspec./m2 g−1 |
---|---|---|---|---|---|
Fe2O3 original powder | 114.1 | 98.8 | 72.9 | 304.9 | 25.0 |
Fe2O3 annealed at 600 °C | 120.0 | 111.7 | 72.2 | 264.7 | 4.9 |
Fe2O3 annealed at 700 °C | 111.8 | 115.9 | 73.6 | 44.8 | 4.4 |
Fe2O3 annealed at 800 °C | 104.3 | 113.9 | 83.9 | 7.9 | 3.0 |
Annealing of the iron oxide powder drastically changes the surface morphology of the material. SEM images of annealed samples are given in Fig. 3. Annealed oxides have a labyrinth structure and become more amorphous from the surface with growth of annealing temperature. Increasing annealing temperature leads to particle growth and to a decrease of the inner particle pores that were in original powder (Fig. 2 and 3). At 800 °C, particles tend to sinter into larger grains and the number of surface pores became minimal (Fig. 3C). Specific surface area of the annealed samples decreases with increasing annealing temperature (4.9, 4.4 and 3.0 m2 g−1; Table 1). Although values of the specific surface area of annealed samples are close to each other, such small differences in morphology affects the TPR spectrum and can be revealed by the analysis of optimized kinetic parameters, as discussed below.
Fig. 3 SEM images of the Fe2O3 samples, annealed at various temperature. (A) 600 °C; (B) 700 °C; (C) 800 °C. |
TPR curves of iron(III) oxide samples, recorded at several heating rates (3, 6, 9 and 12 °C min−1), are presented in Fig. 4A. Using data obtained at several temperature regimes improves the reliability of the results of the kinetic analysis according to the International Congress on Thermal Analysis and Calorimetry (ICTAC) project recommendations.38 Curve fitting was carried out using a three-stage scheme of reduction, with the following kinetic equations:
(20) |
(21) |
Fig. 4 TPR curves of the original Fe2O3 powder (A) recorded under various heating rates, and TPR curves of iron(III) oxide samples (B) recorded under 12 °C min−1 heating rate. |
The shape of the reduction curve changes with heating rate, with process shifts to higher temperatures (Fig. 4A). The first reduction peak shifts by 50 °C from 3 °C min−1 curve to 12 °C min−1. There are also variations in TPR curves if we consider them in the row of studied materials (Fig. 4B). The first maximum in the TPR of iron(III) oxides also shifts to higher temperatures in the progression from original powder to the most annealed. This clearly demonstrates the influence of the surface morphology of the investigated samples on the TPR spectra and, as a result, on the kinetic analysis. However, all TPR curves are approximated well, and the optimization procedure resulted in close values of kinetic parameters demonstrating high reliability of the suggested method (Table 1). It is noteworthy that the preexponential factors found for the first step of reduction correlate well with the specific surface areas of the materials. If we recall that the preexponential factor is a multiplication of the specific surface area at zero conversion and true preexponential factor (9) at a given temperature, then it should correlate to some extent with the specific surface area. Although there is a loose correlation between them, because we deal with the “active” specific surface area that is not the one measured by nitrogen, preexponential factors form the same dependence as the specific surface areas (A1 > A2 > A3 > A4 as well as S1 > S2 > S3 > S4). Thus, we can conclude that the suggested method takes into account the difference in material morphology and is capable of providing reliable kinetic results.
Besides kinetic parameters, we also obtained s(α) curves as presented in Fig. 5. Change in relative specific surface area with degree of conversion for Fe2O3 oxide is a monotone increasing function (Fig. 5A) that can be explained by the contraction particle model. When particles of hematite react with hydrogen, their size decreases and at the same time the specific surface area of the particles increases due to the fact that the ratio between the surface atoms and the bulk atoms grows. The shape of the s(α) curves of the intermediate oxides Fe3O4 and FeO (Fig. 5B and C) is caused by the same factors. The first particles of Fe3O4 formed from the reduced hematite began to grow; thus, their specific surface area decreases. Then, when the reduction rate of Fe3O4 exceeds the growth of Fe3O4, the relative specific surface area begins to grow. It is for this reason that we observe a minimum on the s(α) curve. The same hypothesis we be applied for s(α) dependence of FeO.
Fig. 5 Optimized relative specific surface area vs. conversion degree curves of oxide forms obtained for 12 °C min−1 heating rate. (A) Fe2O3; (B) Fe3O4; (C) FeO. |
From s(α) of Fe2O3, we reproduce the initial particle size distribution using the iterative procedure described in the Calculations section. We evaluated particle size distribution and compared it with the one measured experimentally using a laser diffraction method (Fig. 6). The results have shown that the proposed method allows accurate reproduction of the particle sizes, in particular, the mode in particle size distribution based on very simple theoretical assumptions. Deviations of the model from the experimental values come from the fact that the procedure of nonlinear minimization used in the work is very sensitive to the initial values of the optimized parameters. Good consistency between the model and the theory proves that the suggested method for kinetic analysis of TPR satisfies all the requirements and can be used in wider research practice.
The presented approach was implemented in the kinetic analysis of the TPR of iron(III) oxide samples that possess different surface morphology under different heating rates (3–12 °C min−1). Close values of kinetic parameters were obtained for the studied materials by the presented method, proving its reliability. Correlation between the found preexponential factors and the specific surface areas of the iron(III) oxide samples demonstrates the capability of the method to probe morphology variation among the series of studied materials (from the original Fe2O3 sample to the sample annealed at 800 °C: A* – 304.9 > 264.7 > 44.8 > 7.9; Sspec. – 25.0 > 4.9 > 4.4 > 3.0). Using optimized s(α) dependence, the initial particle size distribution of the sample was obtained and compared with the experimental distribution obtained via laser diffraction. Deviations between model and experimental values can be attributed to the strong dependence of the results from the initial values of the optimized parameters.
Footnote |
† Electronic supplementary information (ESI) available. See DOI: 10.1039/c7ra09848k |
This journal is © The Royal Society of Chemistry 2018 |