The Kinetic Model of Diffusion and Reactions in Powder Catalysts during Temperature Programmed Oxygen Isotopic Exchange Process

The mathematical model of diffusion in powder oxide catalysts during the process of temperature programmed oxygen isotopic exchange is proposed. The diffusion is considered together with the homogeneous and heterogeneous oxygen isotopic exchange processes. The matrix forms of exchange rate equations of simple and complex heteroexchange, and homoexchange reactions which obtain symmetrical forms are analyzed. The quantitative values of model parameters are found from the fitting of experimental data taken from literature of temperature programmed oxygen isotopic exchange process in catalysts ZrO2 and CeO2. The fittings show a good matching of model results with experimental data. The shapes of kinetic curves registered during temperature programmed oxygen isotopic exchange process are analyzed and the influence of various process parameters such as activation energies of simple and complex heteroexchange, oxygen surface concentration of catalyst, ratio of catalysts surface and volume of reactor, diffusion activation energy is considered. The depth profiles of diffusing oxygen species in oxide catalysts powder are calculated.


Introduction
Ceria (CeO 2 ) and zirconia (ZrO 2 ) are the main components of the three-way catalysts used in automobile pollutant abatement [1,2]. They enable a control of the partial pressure of oxygen near the catalyst surface during automotive emission [3,4] and are able to store and/or supply oxygen under fuel-lean and fuel-rich conditions, respectively, what is necessary for the conversion of nitrogen oxides, hydrocarbons and carbon monoxide [5,6]. AO 2 (A = Ce,Zr) releases oxygen during rich conditions and is converted to AO 2−x , whereas during lean conditions AO 2−x is oxidized back to AO 2 [7,8]. These properties are related to the ability of AO 2 system to promote migration/exchange of oxygen species in the reaction [9,10]. Due to high surface and bulk diffusivity, the oxygen atoms can be quickly transferred to the active sites form different parts of the catalysts [11,12]. The highest oxygen storage capacity and catalytic performance of ceria-zirconia materials is at approximately 40-60% of ceria content [13,14].
A powerful tool to study the various oxygen transport processes that may take place in a crystalline oxide is the 18 O- 16 O isotopic exchange technique [15][16][17]. The oxygen isotopic exchange is a suitable method to analyze the diffusivity and interaction of molecular oxygen with a metal oxide. From the registered kinetic curves of oxygen species concentrations, the mechanisms of catalytic processes and heir rates in given catalyst are determined [18,19]. Three types of exchange mechanisms between gaseous dioxygen and lattice oxygen atom have been defined on oxide catalysts [20,21]: (1) homoexchange between adsorbed atoms itself, without participation of atoms form oxide; (2) simple hetero-exchange between one atom of adsorbed oxygen molecule and one atom of the solid oxide and (3) the multiple hetero-exchange when both atoms of adsorbed oxygen molecule are exchanged with solid oxide oxygen atoms. Information on the mechanism of exchange may be obtained from the kinetic curves of oxygen species partial pressures at the beginning of exchange reaction. The formation of 16 O 18 O species as a primary product from 18 O 2 indicates that exchange takes place via a simple hetero-exchange, if 16 O 2 species are formed firstly, it means multiple hetero-exchange mechanism [22]. Isotopic exchange is generally carried out in a recycle, close reactor coupled to a mass spectrometer [23,24]. The evolution of the partial pressures of 16 O 2 , 16 O 18 O and 18 O 2 during the exchange process are registered as a function of time in isothermal conditions. When the range of temperatures where the exchange process occurs is not well-known isothermal conditions could turn out to be tedious. In those cases, the method of temperature-programmed oxygen isotopic exchange is very useful. Using a linear and slow increase of sample temperature, the concentration of the adsorbing/reacting/desorbing molecules are registered [25,26]. The curves obtained in one temperature programing isotopic exchange experiment directly inform on the temperature window where the exchange takes place and give a complete view of the proportion of lattice oxygen atoms (surface and bulk) involved in the process depending of the temperature [27]. Such experiment can permit comparison of different families of oxide samples in an easy way, comparing the evolution of the rate of exchange and the evolution of the number of exchanged atoms versus temperature [28]. Increasing the temperature also enables us to rapidly reach an equilibrium between the 18 O concentration in the gas phase and into the lattice and then to obtain data about the number of exchangeable atoms in the solid [23].
In the temperature programmed isotopic exchange experiments a very important parameter is the surface area of the catalyst. The high surface area provides more active oxygen transfer from catalyst. In this work the influence of surface is of catalyst is investigated. Oxygen transfer form catalysts goes on through the bulk diffusion processes. The bulk diffusion processes in oxide catalysts where was started to investigate in our previous work Ref. [29]. In this work the complete model of bulk diffusion during and exchange reactions in powder catalysts during temperature programmed oxygen isotopic exchange process is proposed and verification with experimental results is achieved.

The Model
Isotopic oxygen exchange on surfaces of catalysts goes on by homoexchange and by heteroexchange [30,31]. Homoexchange occurs without participation of surface oxygen. The heteroexchange takes place when oxygen atoms form the surface of oxide are involved into the process. The simple heteroexchange occurs when one oxygen atom in the molecule is replaced, in complex heteroexchange both oxygen atoms in molecule are replaced [18,22,32]. Simple and complex oxygen isotopic heteroexchange reactions are listed in Table 1. Not all of the simple hetero-exchange reactions listed in Table 1 result in a change in system composition. It may happen that an atom of the same type from the gas phase can be replaced with an atom of the same type from the surface (e.g., 16 O g with 16 O s ). In that case microscopically the oxygen exchange occurs, but the composition of system remains the same. In order to involve into the calculation only those cases of reaction which lead to compositional changes of system, in Table 1 for simple heteroexchange reactions the probabilities are indicated. Through the probabilities that are included only those cases of exchange when composition of system occurs. Considering adsorption, in the proposed model it is assumed that adsorption is very fast process [33,34] and surface is covered by oxygen species with the same composition as composition in gas phase. The change in composition in the gas phase and on the surface is calculated according to the law of mass action well known in chemistry. Using this law in the case of simple heteroexchange the variation of concentrations in gas phase n 32 , n 34 and n 36 of species 16 where V is volume of reactor and S is surface area of catalyst. Only equation for c 32 is symmetrical to equation for n 34 . Ratio V/S appears because of different dimensions of concentrations in gas phase n ij expressed in mol/m 3 and surface concentrations c ij expressed in mol/m 2 . Later it will be shown that this ratio, which seems to be just technical parameter in fact plays very important role. Using law of mass action in the case of complex heteroexchange the variation of concentrations in gas phase is written as follows: where k C is reaction rate constant of complex heteroexchange expressed by Arrhenius law k C = A C exp(−Q C /RT), where A C and Q C are pre-exponential term and activation energy of complex heteroexchange. Due to the symmetry of complex heteroexchange reactions (that is not the case for simple heteroexchange), the variations of variables c 32 , c 34 and c 36 can be written in very simple form.
From the mathematical point of view, it is very interesting to analyze obtained sets of Equations (1)- (3). For example, the Equation (3) keeping n 32 , n 34 and n 35 as variables can be written in the following matrix.
The Equation (5) is symmetrical in various cross sections and show interesting mathematical regularities of system where the complex heteroexchange of isotopes takes place. For simple heteroexchange the Equation (1) obtains the following matrix form.
Equation (6) also is symmetric in various cross sections and show mathematical regularities in complex heteroexchange systems. More detail and deep analysis of those matrixes could give very useful additional information about the properties of the system.
In the case of homoexchange, when exchange occurs only between adsorbed molecules without participation of lattice atoms the reactions in which the change of system composition takes place are following [20,21].
The rate equations mathematically describing those reactions applying the mass action law are the next:  (5) and (6), if writing Equation (8) in matrix form keeping n 32 , n 34 and n 35 as variables the following two equivalent matrixes can be written Considering the process of bulk diffusion, when oxygen isotope atoms 18 O penetrate in deeper layers of oxide, the concentration variation of oxygen atoms in one oxide layer K as given is calculated by using the second Fick's law expressed in finite increments. However, first it must be adapted to geometrical specifics of particles of powder catalysts. The particles of ceria and zirconia powder catalysts, which experimental results will be fitted by proposed model, have cubic-like geometrical form [35,36]. Describing mathematically the process of diffusion into bulk of powder particle it is necessary to assume the limitation of the depth. The Fick's law expressed in finite increments means that the bulk is deleted into separate layers but in the cases of powder particles, when diffusion flux takes place from for all surfaces into the center of powder particle the area of layers decreases. To solve this problem, the cubic-like powder particles are virtually divided into four pyra-mids and in 2-d case each of them is divided into layers (see Figure 1). The area of each K layer can be found form the following relation: where S ox is the area of surface layer, h is the thickness of layer and d ox is the size of powder particle. The second Fick's law expressed in finite increments, assuming decreasing areas of layers obtains the following form: where D is the bulk diffusion coefficient of oxygen atoms in oxide, i = 16, 18 indicates the type of oxygen isotopes and coefficients B (K) involve the changes of areas of layers and are expressed as: It is assumed that diffusion of oxygen isotopes are balanced, i.e., 18 O diffuses into the bulk and replaces 16 O atoms which diffuses to the surface and the condition is fulfilled. Diffusion coefficient of both oxygen species is assumed the same. As indicated in Equation (11) K > 1 the equation describes diffusion starting from the second layer K ≥ 2, K = 1 is the surface, first or adsorption layer where oxygen species form gas phase adsorbs (and desorbs). The equation is needed which could describe diffusion and oxygen exchange between surface (first) layer and second layer of catalysts. The mass action law in combination with Fick's law we will use to build equation for diffusion between first and second monolayers.
The exchange reactions between oxygen species on the surface 18 Using the mass action law and taking into account the gradient of concentrations the reaction rates of those reaction are expressed in following form: where k d is the reaction rate constant, c 16 and c Finally, the kinetics of composition in gas phase is calculated including all considered above processes simple and complex heteroexchange, homoexchange and diffusion: The temperature programming exchange.
In the model the reaction rate constants k S , k C , k O , k d and diffusion coefficient D depend on temperature according to Arrhenius law k j = A j exp −Q j /RT(t) , A O and Q O are preexponential term and activation energy of. In temperature programming exchange process the temperature T depends on time linearly: T(t) = T 0 + bt, where T 0 is initial temperature and b the rate temperature increase, t is the time.

Results and Discussion
In Figure 2a,b the experimental and calculated dependencies of partial pressures of oxygen species 18 [35,36]. The temperature programmed oxygen isotopic exchange experiment was performed at following conditions: the temperature increase speed 1.6 • C/min, specific surface area S bet = 25 m 2 /g for ceria and S bet = 25 m 2 /g for zirconia, mass of powder catalyst m cat = 0.25 g, initial pressure of 18 O 2 48 mbar and volume of reactor V = 12 cm 3 . The surface concentration of oxygen in ceria depending on crystallite orientation varies from 0.97 nm 2 for (110) until 15.8 nm 2 for (111) [35,36] and in these calculations was taken as 1.09 × 10 19 m 2 . Similarly using data from Ref. [36] for zirconia oxygen surface concentration was taken 1.09 × 10 19 m −2 . From the fitting of experimental curves presented in Figure 2, the values of activation energies of simple and complex heteroexchange reactions, diffusion and its pre-exponential terms were found. Values are written in Table 2.
The fitting results presented in Figure 2 are quite good, taking into account the difficulty of the process. In TPIE experiments the temperature changes with time and curves in Figure 2 represent not a simple dependency on temperature but the kinetic curves. Some deviations of calculated results form experimental ones in many cases occurs because of impurities, which always exist in real conditions and which are impossible to estimate in modeling. For zirconia Figure 2a some deviation of calculated curves form experimental points occurs considering the temperature range where the process of exchange starts. Experimental results show that exchange process starts at a little lower temperatures tha, theoretical predictions. However, for ceria catalysts Figure 2b the calculated curves very well correspond with experimental points in this interval of temperatures. For ceria a small deviation is observed at higher temperatures when steady state of process is reached, but only for p 34 curve, the curve p 32 show very good agreement in whole considered temperature interval. Comparing the temperature interval where the exchange process starts the lower temperature is for ceria catalyst.   Figure 2.

Parameters Ceria Zirconia
Measuring composition changes in gas phase during TPIE process the important parameter is the ratio between volume of reactor and total surface of catalyst V/S which is involved into calculations through Equations (2) and (4). This ratio influences significantly the partial pressures of oxygen species, including the steady state regime. This ratio depends on mass of catalysts and specific surface area through the relation S = m cat ·S bet . This relation also means that the same influence in the kinetic curves have both parameters mass of catalyst and specific area of catalyst. Moreover, because in Equations (2) and (4) the ratio V/S exist, the volume of catalyst also significantly affects the shapes kinetic curves and the values of partial pressures when steady state is reached. At the steady state if when mass of catalysts is relatively high p 34 < p 32 and when catalyst mass is low, situation is in opposite: p 34 > p 32 . However, it is important to consider the influence of process parameters on the shapes of kinetic curves and partial pressures of species at steady state regime.
In Figure 3 the calculated curves of partial pressures of oxygen species at steady state regime as a dependency on ratio S/V are presented. The shapes of curves p 34   Experimental kinetic curves of oxygen species partial pressures registered during TPIE may have quite different shapes. For example, the curve p 34 in Figure 2b has a maximum, but the curve p 34 not. Many experimental results show that the maximums can be broad or narrow, high or low, or sometimes they are not formed at all. It depends on type of catalyst. In order to clarify this situation and to find some regularities the calculations were performed by varying activation energies of simple and complex heteroexchange. The influence of Q S and Q C (activation energies of simple and complex exchange) is analyzed in Figure 4a the calculated curves obtained by changing Q C and (b) the curves by changing Q S are presented. From Figure 4a it is seen that with increase of Q C the temperature at which the molecules 16 O 2 start to form shifts to higher temperatures and this shift at the beginning is very significant and later reduces (see 4 and 5 curves). However, the temperature of formation of 16 O 18 O species almost does not depend on Q C but the shape of p 34 curves significantly depends on Q C . At low values of Q C the maximum is not observed in those curves, but it appears and increases with the increase of Q C . In opposite it is for p 32 curves, at low Q C the broad but not well-expressed maximum is seen but with increase of Q C it disappears. The influence of Q S which presented in Figure 4b is different. The temperature of formation of 16 Figure 4a,b it is seen that Q S and Q C values determine which oxygen species with increase of temperature will be formed first. Depending on Q S and Q C values with increase of temperature the oxygen species 16 O 18 O can start to form first and 16 O 2 after or in opposite, 16 O 2 first and 16 O 18 O after. The cases when both oxygen species start to form together at the same temperature also can be observed in Figure 4. The results presented in Figure 4 were obtained at relatively high ratio of S/V, in the case when steady state pressures as p 32 > p 34 . The case of relatively low ratio S/V when steady state pressures are p 32 < p 34 is considered in Figure 5. In Figure 5a the influence of Q C is shown. It is seen that curves of p 34 almost are not affected by changes of Q C , but curves p 32 are significantly affected by Q C . At low values of Q C curves p 32 have big and broad maximums which decrease and finally disappear with increase of Q C . In Figure 5a it is also seen that at low values of Q C with increase of temperature species of 16 Figure 5b. Now, the change of Q S significantly affects both curves. At higher values of Q S the maximums in curves p 32 start to appear and become broader with increase of Q S . However, it is interesting to note that when maximum appears, the temperature of formation of 16 O 2 does not change any more with further increase of Q S . The curves of p 34 continuously shift to higher temperatures with increase of Q S but the shape of curves remain the same and no maximums appear. At low values of Q S with increase of temperature the species of 16   Considering the steady state pressures of oxygen species, it is important to note that additionally to ratio S/V, there is another parameter which influences the steady state pressures and modifies the kinetic. This parameter is the initial surface concentration of oxygen C OX . In the model this parameter is involved through the initial value of surface concentration c 32 . Parameter c 32 is variable but at the beginning at t = 0 it equals to C OX : c 32 (t = 0) = C OX . Surface concentration of oxygen depends on the type of catalyst, method of preparation and also orientation of oxide crystallite grains [35,36]. The influence of oxygen surface concentration is shown in Figure 6 where the partial pressure curves calculated for different surface concentration of oxygen C OX are presented. It is seen that in contrast with curves presented in Figures 4 and 5, the curves of Figure 6 at steady state regime differs. Since C OX depends on crystallite surface orientation curves presented in Figure 6 can also be considered as presentation of the influence of crystallite surface orientation catalysts. Parameter C OX affects both, the temperature of formation of oxygen species and the pressure at steady state regime. In the heterogeneous exchange, the process of diffusion is very important and influences the kinetic curves. The influence of bulk diffusion process is shown in Figure 7 where the calculated kinetic curves of partial pressure of oxygen specie 16 O 2 as a dependency on temperature (which increases with time in TPIE) are presented. The curves are calculated at different values of diffusion activation energy at constant pre-exponential factor A dif = 1.90·10 −24 m 2 s −1 . All other parameters were kept as constants and were the same as in calculations presented in in Figure 2b. It is seen that diffusion significantly influences the shape of kinetic curves in the transition period and at steady state. The partial pressure at steady state of 16 O 2 decreases with increase of Q dif . This result is logic, because the increase of Q dif means the decrease of diffusivity. At higher diffusivity more amount of 16 O 2 appear in gas phase because of diffusion of oxygen atoms from the bulk of catalyst. In order to see this effect in more detail the dependence of 16 O 2 partial pressure at steady state regime is presented in Figure 8. The considered interval is narrow and change of pressure is small, nevertheless, the nonlinear dependence of steady state partial pressure on diffusivity can be seen.
In order to check the validity of model the concentration distribution in deeper layers, the concentration deps profiles of atoms 16 O and 18 O were calculated. Results are presented in Figure 9, where depth profile curves are calculated for different values of Q dif . Obtained curves are typical diffusion curves and shows the correctness of calculations. However, the bulk diffusion on powder catalysts is not so simple as discussed above presenting Figure 1 and Equation (12) and depth profile curves at higher diffusivity can be more complex. It is seen that with increase of diffusivity (decrease of Q dif ) the concentration of 16 O at the surface layers decreases and is replaced by oxygen 18 O coming from gas phase. It means what more oxygen 16 O must appear in gas phase forming oxygen species 16     Above the bulk diffusion was considered, but the surface diffusion is also very important in catalysis. It depends on homogeneity of surface. In the case of homogeneous surface, the exchange takes place directly on surface of oxide, the surface diffusion does not change surface composition of species and mathematically there is no need to write new equation because of absence of new variables. Concentration gradients on the surface in that case are not formed and it is not possible to consider surface mobility of atoms. Situation is significantly different in the case of nonhomogeneous surfaces, e.g., when noble metal nanoparticles are formed on the surfaces of oxides such as in M/CeO 2 (M-noble metal, e.g., Pt, Pd, Au) catalysts. Such catalysts are used in order to reduce oxygen exchange temperature [23] (which on noble metals is less) what is very important in automotive catalysts. In that case, because of spillover [12,37,38] the oxygen concentration gradients are formed toward metal nanoparticles. The model of surface and bulk diffusion because of oxygen spillover was proposed in our previous works [33,34]. In these studies the surface diffusion was combined with bulk diffusion and finally two dimensional diffusion model was created for consideration of nonhomogeneous catalysts such as M/Ce x Zr (1−x) O 2 . It was found [28,39] that process of bulk diffusion becomes more important with increase of content of Zr. For pure ceria CeO 2 samples bulk diffusion is weak and dominates surface diffusion.

1.
The matrix forms of equations of simple and complex isotopic heteroexchange have symmetry in various cross sections and show interesting mathematical regularities.

2.
Despite the variety of different shapes of partial pressure kinetic curves registered during TPIE there are strict regularities: (a) in the case of p 34 < p 32 at steady state, the p 34 curves have well expressed maximums or may be without extremums at higher values of activation energies or together with maximum can exist minimum; the p 32 curves approaches steady state passing broad maximums or without them but never through the minimums. (b) in the case p 34 > p 32 at steady state, the p 32 curve approaches steady state passing broad maximums or without them at higher activation energies, p 34 curve approaches steady state without extremums.

3.
The main parameter determining the steady state partial pressures of oxygen species is the ratio of catalyst specific surface area to volume of reactor. Surface atomic concentration bulk diffusion also influences state partial pressures. Steady state partial pressures are not influenced by reaction rates of heteroexchange.

4.
The bulk diffusion process significantly influences the shapes of kinetic curves and partial pressures of species at steady state regime.