Literature Review on Thermodynamic and Kinetic Limitations of Thermal Decomposition of Methane

Abstract

The state of the art in methane pyrolysis does not yet provide a definitive answer as to whether the concept of an elementary reaction is universally applicable to the apparently simple process of methane dissociation. Similarly, the literature currently lacks a comprehensive and unambiguous description of the methane pyrolysis process and, in particular, a single model that would well represent its course at both the micro and macro scales. Given the wide range of conditions under which this reaction can occur—whether thermal or thermo-catalytic, in solid or fluidized bed reactors—it is crucial to evaluate the usefulness of different kinetic models and their compatibility with basic thermodynamic principles and design assumptions. To address these research gaps, the authors analysed the thermodynamic and kinetic dependencies involved in the thermal decomposition of methane, using the synthesis of methane from its elemental components and its reversibility as a basis for exploring suitable kinetic models. Using experimental data available in the literature, a wide range of kinetic models have been analysed to determine how they all relate to the reaction rate constant. It was found that regardless of whether the process is catalytic or purely thermal, for temperatures above 900 °C the reversibility of the reaction has a negligible effect on the hydrogen yield. This work shows how the determined kinetic parameters are consistent with the Kinetic Compensation Effect (KCE) and, by incorporating elements of Transition State Theory (TST), the possibility of the existence of Entropy–Enthalpy Compensation (EEC). The indicated correspondence between KCE and EEC is strengthened by the calculated average activation entropy at isokinetic temperature $({∆S}_B=-275.0 \mathrm{J}·(\mathrm{m}\mathrm{o}\mathrm{l}·\mathrm{K}{)}^{-1}).$ Based on these results, the authors also show that changes in the activation energy $(E=20–421 \mathrm{k}\mathrm{J}·{\mathrm{m}\mathrm{o}\mathrm{l}}^{-1})$ can only serve as an estimate of the optimal process conditions, since the isoconversion temperature ($T_{iso}=1200-1450 \mathrm{K}>T_{eq}$) is shown to depend not only on thermodynamic principles but also on the way the reaction is carried out, with temperature (T) and pressure (P) locally compensating each other.

1. Introduction

Among the known C₁ to C₁₄ n-aliphatic hydrocarbons, methane is one of the most stable (${∆}_{298} G^{\varnothing }=-50.76 \mathrm{k}\mathrm{J}·{\mathrm{m}\mathrm{o}\mathrm{l}}^{-1}$) [1]. However, the introduction of oxygen into the reaction system renders methane one of the most reactive of the aliphatic hydrocarbons. It is noteworthy that methane, when reacting with the source of oxygen (any form of oxygen, e.g., oxides), can lead to the production of not only CO₂ and CO, but also coupling products such as C₂ hydrocarbons (ethane, ethene) and oxygenated derivatives (formaldehyde, methanol) [2].

The current shift in the structure of energy consumption is, among other things, leading to an intensification of the search for hydrogen production methods. This situation has led to an interest in the process of methane pyrolysis. In principle, the aim of this process is comparable to the well-known steam–methane reforming, i.e., the production of hydrogen. The oxygen present in the environment of the steam–methane reforming process greatly enhances the aforementioned reactivity, but the former has the theoretical advantage of being inherently free of CO₂ emissions. Methane pyrolysis is an endothermic process (${∆}_r H^{\varnothing }=74.85 \mathrm{k}\mathrm{J}·{\mathrm{m}\mathrm{o}\mathrm{l}}^{-1}$) [3] which occurs due to the decomposition of the methane molecule into its constituent elements, solid carbon and gaseous hydrogen. The following notation is widely accepted as the reaction of the decomposition of methane:

$$\mathrm{C}{\mathrm{H}}_4\rightleftarrows \mathrm{C}+2{\mathrm{H}}_2$$

(1)

The efficiency of reaction (1) can be influenced by changing the process parameters such as temperature, pressure, and the introduction of a catalyst. Another key factor determining the course of reaction (1) is the method of heat supply. A general distinction is made between indirect and direct heat supply, but the type of reactor used can also have a major influence on the process. From the point of view of reactor design, it is known that fixed bed, fluidised bed, or liquid bed reactors are the most commonly used. In addition, when considering temperature ranges > 1300 °C or processes where plasma is used as a heat source for the reaction, the reaction mechanism (1) may change and other products besides hydrogen and structural carbon materials can be expected [3,4,5]. If reaction (1) is conducted at temperatures corresponding to the plasma state, upon quenching the resulting gaseous products, the following reactions may take place:

$$2\mathrm{C}{\mathrm{H}}_4\rightleftarrows \mathrm{C}\mathrm{H}\equiv \mathrm{C}\mathrm{H}+3{\mathrm{H}}_2$$

(2)

$$4\mathrm{C}{\mathrm{H}}_4\rightleftarrows \mathrm{C}\mathrm{H}\equiv \mathrm{C}\mathrm{H}+\mathrm{C}{\mathrm{H}}_2=\mathrm{C}{\mathrm{H}}_2+5{\mathrm{H}}_2$$

(3)

Thus, reaction (1) should only be considered as a conceptual/general notation of the thermal decomposition of methane, while reactions (2) and (3) can illustrate examples of pathways for the formation of acetylene from methane. It should be noted that acetylene can also be obtained from methane by partial combustion, while the trimerization of acetylene is a pathway leading to the synthesis of benzene and higher aromatic hydrocarbons. A graphical representation of the network of reactions that can lead from methane to C_5–6 hydrocarbons as possible products of its pyrolysis is given in [6].

As mentioned above, the course and efficiency of reaction (1) can be influenced by many factors. One of the most important is the characteristics of the catalyst used. In [7] most of the metals known for their activity in methane decomposition were ranked in decreasing order of their influence on the efficiency/rate of reaction (1): Ni, Co, Ru, Rh > Pt, Re, Ir > Pd, Cu, W, Fe, Mo. Due to their abundance, cost, environmental and process characteristics, the greatest attention is currently being paid to the catalysts based on Ni, Co and Fe. Belonging to the group of transition metals, these catalytic materials have partially filled 3d orbital [7], which is an important similarity with the catalysts used in industrial Fischer–Tropsch synthesis. It is precisely because of the possibility of tailoring the properties of the $\mathrm{C}$ product that the literature on the selection of catalysts for reaction (1) is so extensive. Among the above-mentioned catalysts, nickel [3,4,5,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,23,24,25,26,27,28,29,30,31,32,33,34,35,36,37,38,39,40,41,42,43,44,45,46,47] occupies a special position due to its effectiveness and the possibility of tailoring the C-structures, while iron [3,4,5,7,8,27,31,40,42,48,49,50,51,52,53,54,55,56,57,58,59,60,61] is notable for being much cheaper and more environmentally friendly.

It is worth noting that iron ore, which can be used as a catalyst for methane decomposition, has a unit cost of about 100 times lower than nickel catalysts [57]. Another argument in favour of iron catalysts is the possibility of using the spent catalyst, for example in metallurgical processes (blast furnace, foundry), which, when decomposing biomethane, can lead to an effective method of decarbonising this ‘hard to abate’ sector. As a last resort, the spent iron catalyst can be easily and safely disposed of, which is not the case with Ni. The main disadvantage of Fe compared to Ni is its lower activity and therefore the need to use higher temperatures. While theoretically 600–800 °C is sufficient to make $\alpha -Fe$ catalytically active [49] even if the process were carried out under fixed bed conditions and at 800 °C, it would still be technically infeasible due to the time required for decomposition [57]. Importantly, Fe also has a very sharp upper-temperature limit due to the polymorphic transformation of iron that occurs at 912 °C.

The characteristic that unites both iron and nickel in the context of reaction (1) is the multifunctional catalytic role of these materials, and it is noteworthy that this characteristic is commonly attributed to the formation and decomposition of carbides. In this respect, the processes associated with the carburisation and decarburisation of iron [62] deserve particular attention, which is the reason why so much work has been performed to date on the formation/decomposition reactions of cementite (Fe₃C) [60,63,64,65,66,67] and Hägg carbide (Fe₅C₂) [68]. In the context of iron, the reactivity and properties of cementite are considered to be of significant importance for its catalytic activity [68,69,70,71,72,73], while for nickel similar materials can be found for nickel carbides such as Ni₃C, Ni₈C [13,45], or Ni₈C on TiC substrate [47]. Importantly, the literature for the latter is much sparser.

The above series did not cover the last material of importance in terms of catalysis of reaction (1), i.e., internally structured carbon materials. The group of carbon catalysts finds somewhat less attractiveness in methane pyrolysis, mainly due to their lower catalytic activity, the requirement to perform the process at temperatures ≥ 1000 °C (for reference 750 °C for Ni and 850 °C for Fe) and the lower possibility to adjust the quality of the produced carbon structures.

All available data agree that regardless of the type of catalyst used, the deposition of C-structures on its surface is the main reason for deactivation, but since some forms of carbon deposits show some limited catalytic activity, the process reaction (1) can theoretically become autocatalytic for a very specific range of parameters. A drawback of this reasoning may be the vague distinction between a purely thermal or C-catalysed process at temperatures > 1200 °C.

Depending on the process conditions, including the characteristics of the catalyst used, the carbon structures produced in reaction (1) consist of either a dominant or a mixture of different structural conformations. Thus, the literature describes carbon black, mesoporous carbons, carbon nanotubes (CNTs), carbon nanofibers (CNFs), multi-walled carbon nanotubes (MWCNTs), amorphous turbostratic carbon, graphite-like carbon, filamentous (or foam) carbon as solid products of the methane pyrolysis reaction. Together these materials are referred to as the ‘carbon family’.

When considering the course of a catalysed reaction (1), the emergence of a new solid phase interacting with the process induces a much higher level of complexity in the analysis. According to one of the most frequently discussed mechanisms [4], molecular or dissociative adsorption on catalysts involves the formation of radicals: CH₃^•, CH₂^•, CH^•, C^•, H^• and their subsequent reduction and recombination. The final steps in the molecular or dissociative adsorption models are the nucleation of free carbon and the growth of its deposit [5]. It is generally accepted that the radical mechanism dominates the reaction carried out at temperatures > 1200 °C [41]. Further details on reaction mechanisms can be found in [7,12,65,74,75,76,77,78,79].

2. Aim of the Work

The state of the art in methane pyrolysis does not yet provide a definitive answer as to whether the concept of an elementary reaction is universally applicable to the apparently simple process of methane dissociation. Similarly, a comprehensive and unambiguous description of the methane pyrolysis process is currently lacking in the literature and, in particular, a single model that would well represent its course at both micro and macro scales. In order to address these issues, the authors present below an analysis aimed at determining whether it is possible to fulfil this purpose using an approach typical of phenomenological thermodynamics. Therefore, in the following parts of this work, the applicability of using kinetic equations as a tool to describe the course of reaction (1) was verified. In this way, the correctness of using the concept of an elementary reaction in the context of the methane pyrolysis process will also be addressed. In order to test the hypothesis that the constraints (limits) resulting from phenomenological thermodynamics alone are sufficient to verify the above considerations, the authors have taken advantage of experimental results presented in the already published literature. It is important to note that the Kinetic Compensation Effect (KCE) and the equivalent elements of the Transition State Theory (TST), expressed by the Enthalpy–Entropy Compensation (EEC), are of paramount importance for this approach.

3. Phenomenological Thermodynamics and Thermodynamics of Activation

From a thermodynamic point of view reaction (1) is reversible and is expressed by the equilibrium constant, which is a function of temperature and pressure. Its equation can be written as [3,4]:

$$K=K_p{\left({\displaystyle \frac{P^{\mathrm{\varnothing }}}{P}}\right)}^{∆\nu }$$

(4)

This equation is further simplified for the case of reaction (1) and the gaseous components, since then $∆\nu =1$.

The number of terms in the polynomial and the values of its coefficients with respect to temperature depend on the assumed range of variation and the expected accuracy, assuming a gas pressure activity of 1. An extended expression is obtained from the free enthalpy versus temperature equation presented in the paper [20]:

$${∆}_r{G}^{\varnothing }=89.66-0.1023·T-4.28·{10}^{-6}{T}^2-{\displaystyle \frac{2499.36}{T}},\mathrm{k}\mathrm{J}·{\mathrm{m}\mathrm{o}\mathrm{l}}^{-1}$$

(5)

By dividing both sides of Equation (5) by $\left(-RT\right)$, we obtain:

$$-{\displaystyle \frac{{∆}_r{G}^{\mathrm{\varnothing }}}{RT}}=-{\displaystyle \frac{89.66}{RT}}+{\displaystyle \frac{0.1023}{R}}+{\displaystyle \frac{4.28·{10}^{-6}T}{R}}+{\displaystyle \frac{2499.36}{R{T}^2}},$$

and given that $-{\displaystyle \frac{{∆}_r{G}^{\mathrm{\varnothing }}}{\mathit{RT}}}=\ln K_p$, this expression can be further simplified to the form:

$$\ln K_p=-{\displaystyle \frac{10784.22}{T}}+12.30+0.000515·T+{\displaystyle \frac{300620.64}{T^2}}$$

(6)

Solving Equation (6) for $∆G^{\varnothing }=0$ and ${ K}_p=1$, we can determine that $T_{eq.}=818.84\mathrm{K}$.

A modified, shorter version of Equation (5) can be expressed using the relationship given in [55].

$${∆}_r{G}^{\mathrm{\varnothing }}=88.039-0.1079·T,\mathrm{k}\mathrm{J}·{\mathrm{m}\mathrm{o}\mathrm{l}}^{-1}$$

(7)

In a similar way to the recalculation presented above for Equation (5), if both sides of Equation (7) are divided by $\left(-RT\right)$, then we obtain

$$\ln K_p=-{\displaystyle \frac{10589.25}{T}}+12.978$$

(8)

Equation (8) is applicable to similar temperature and pressure ranges (${573 \mathrm{K}\le T}_{eq}\le 1773 \mathrm{K}$, $P=0.1 \mathrm{M}\mathrm{P}\mathrm{a}$), by solving it for $∆G^{\varnothing }=0$ and ${ K}_p=1$, we can determine that $T_{eq}=815.93 \mathrm{K}$. This simple exercise shows us that, for reaction (1), we can obtain comparable results even by using very simple approaches determined from different angles of analysis. Considering the constants mentioned above, a very similar relationship can be derived from [80]. For this reason, Equations (4)–(8), and in particular Equations (7) and (8), are used as the basis for the following discussion.

On the other hand, algebraic forms of thermodynamic equilibrium constants ${K=K}_p{=K}_x{=K}_{\alpha } \mathrm{e}\mathrm{t}\mathrm{c}.$ (for constant total pressure) are correct for the following experimental (without subscript) and equilibrium (with the index eq) values.

$${K=K}_x={\left.{\displaystyle \frac{{x_{{\mathrm{H}}_2}}^2}{1-x_{{\mathrm{H}}_2}}}\right|}_{eq}={\left.{\displaystyle \frac{(1-x_{\mathrm{C}{\mathrm{H}}_4}{)}^2}{x_{\mathrm{C}{\mathrm{H}}_4}}}\right|}_{eq}$$

(9)

Furthermore because

$$x_{{\mathrm{H}}_2}+x_{\mathrm{C}{\mathrm{H}}_4}\equiv 1$$

(10)

for the condition of constant flow of methane $\dot{V}>0$, the conversion of methane can be written as follows:

$${\alpha }_{{\mathrm{C}\mathrm{H}}_4}={\displaystyle \frac{x_{{\mathrm{H}}_2}}{2-x_{{\mathrm{H}}_2}}}={\displaystyle \frac{1-x_{\mathrm{C}{\mathrm{H}}_4}}{1+x_{\mathrm{C}{\mathrm{H}}_4}}}$$

(11)

Using Equations (10) and (11) the thermodynamic equilibrium constants can be expressed also in terms of the conversion degree of methane:

$${K=K}_{\alpha }={\left.{\displaystyle \frac{4{{\alpha }_{{\mathrm{C}\mathrm{H}}_4}}^2}{1-{{\alpha }_{{\mathrm{C}\mathrm{H}}_4}}^2}}\right|}_{eq}$$

(12)

Similarly, relations (10) and (11) are correct for both experimental and equilibrium data derived from thermodynamic equations. More detailed thermodynamic considerations of reaction (1), including its mechanism, can be found in the following works [11,13,42,65,74,78,81,82,83,84].

Phenomenological thermodynamics can only provide useful information in terms of the constraints imposed by the equilibrium states of a given reaction. Since the course of reaction (1) is facilitated by the use of catalysts, it must be analysed in the context of the reactions and transformations taking place on the surface of the catalyst, so that the interactions between the reactants and the surface morphology and chemistry of the catalyst are also of great importance. For the reasons given in the introduction above, the following considerations are based on iron as the catalyst with probably the best all-round properties.

In the case of iron catalysts, particular attention should be paid to iron ores, as these materials can provide a good balance between the optimum process temperature, low purchase and preparation costs and the quality of the carbon structures produced. From the family of ferrous materials, the following minerals are of particular interest: hematite, magnetite, goethite, and siderite. It is noteworthy that, despite physical differences, their activity can typically be described by comparison with the characteristics of haematite ($\mathrm{F}{\mathrm{e}}_2{\mathrm{O}}_3$) or $\mathsf{\alpha }-\mathrm{F}\mathrm{e}$ [85].

In a highly reducing atmosphere, a series of reduction reactions are observed. Their course and extent depend on the temperature and the chemical properties of the reducing agent [40,73,85,86,87,88]:

$$\mathrm{F}{\mathrm{e}}_2{\mathrm{O}}_3\to \mathrm{F}{\mathrm{e}}_2{\mathrm{O}}_3·\mathrm{F}\mathrm{e}\mathrm{O}\to \mathrm{F}{\mathrm{e}}_{\mathrm{x}}\mathrm{O}\underset{>575\mathrm{℃}}{\to }\mathrm{F}\mathrm{e}\mathrm{O}\underset{<912\mathrm{℃}}{\to }\mathsf{\alpha }-\mathrm{F}\mathrm{e},$$

(13)

where $x=0.82÷0.95$. When the atmosphere is less reducing, such as in the case of dilute methane, haematite is reduced to wüstite ($\mathrm{F}{\mathrm{e}}_{\mathrm{x}}\mathrm{O}$) at a much faster rate than subsequently to $\mathrm{F}\mathrm{e}$ [87], and this process occurs by two parallel mechanisms [88]. Similar phenomena have also been observed in the case of hydrogen reduction, as the reduction in FeO is very slow [40]. It is noteworthy that further comparative studies on the methane reactivity of the above-mentioned pair of Fe oxides are supported mainly by the different reaction rates [41,89].

Allotropic forms of iron are stable within certain temperature ranges [85]:

$$\mathsf{\alpha }-\mathrm{F}\mathrm{e}<912 °\mathrm{C}, \mathsf{\gamma }-\mathrm{F}\mathrm{e} 912–1394 °\mathrm{C}, \mathsf{\alpha }\left(\mathsf{\delta }\right)-\mathrm{F}\mathrm{e} 1394–1538 °\mathrm{C}$$

(14)

where $\mathsf{\alpha }-\mathrm{F}\mathrm{e}$ has a body-centred cubic (BCC) structure, $\mathsf{\delta }$ represents the transformation of $\mathsf{\alpha }-\mathrm{F}\mathrm{e}$ to a larger lattice size, and $\mathsf{\gamma }-\mathrm{F}\mathrm{e}$ has a face-centred cubic (FCC) structure.

A comparison of the heat capacity for the polymorphic transformation $\mathsf{\alpha }-\mathrm{F}\mathrm{e}$ → $\mathsf{\gamma }-\mathrm{F}\mathrm{e}$ is shown in Figure 1. While the maximum of the blue curve indicates the Curie point ($T=1042 \mathrm{K})$, after integration and recalculation of this data to the unit of $\mathrm{k}\mathrm{J}·\mathrm{m}\mathrm{o}{\mathrm{l}}^{-1}$, the following relationship for the Gibbs free energy of $\mathsf{\alpha }-\mathrm{F}\mathrm{e}$ can be determined for the temperature range of $700–1300 \mathrm{K}: G^{\varnothing }=-0.070·T+19.65$. The thermal characteristic of $\mathsf{\gamma }-\mathrm{F}\mathrm{e}$ can be further completed by the fact that the linear characteristic presented here disappears at the melting temperature of iron. Finally, the polymorphic transformation of $\mathsf{\alpha }-\mathrm{F}\mathrm{e}$ is endothermic and takes place in the temperature range $1000–1200 \mathrm{K}$. The enthalpy of this transformation has the value of $∆H^{\varnothing }=10.29 \mathrm{k}\mathrm{J}·\mathrm{m}\mathrm{o}{\mathrm{l}}^{-1}=184.24 \mathrm{k}\mathrm{J}·\mathrm{k}{\mathrm{g}}^{-1}$, while for reference, the standard enthalpy for the decomposition reaction of methane is $∆H^{\varnothing }=74.85 \mathrm{k}\mathrm{J}·\mathrm{m}\mathrm{o}{\mathrm{l}}^{-1}$.

Table 1. Compilation of the chemical reactions and their subsequent parameters of the Gibbs free energy equations [${∆}_{r}G=a+bT$, $\mathrm{k}\mathrm{J}·\mathrm{m}\mathrm{o}{\mathrm{l}}^{-1}$ $P=0.1 \mathrm{M}\mathrm{P}\mathrm{a}$] of the species present in the atmosphere of the methane pyrolysis process when iron oxides are used as catalysts.

No.	Reaction	a, kJ·mol−1	b, kJ·(mol·K)−1	T Range, K	${\mathit{T}}_{\mathit{e}\mathit{q}}$$\mathbf{for} {∆}_{\mathit{r}}\mathit{G}=0$, K	Source	Remarks
1	$\mathrm{C}{\mathrm{H}}_4\rightleftarrows \mathrm{C}+2{\mathrm{H}}_2$	88.04	−0.108	573–1773	815.2	[55]	Equation (7)
2	$2{\mathrm{F}\mathrm{e}}_2{\mathrm{O}}_3+3\mathrm{C}{\mathrm{H}}_4\rightleftarrows$ $6{\mathrm{H}}_2+3{\mathrm{C}\mathrm{O}}_2+4\mathrm{F}\mathrm{e}$	248.37	−0.289	298–1300	859.4	[1] *	[59]
3	${\mathrm{F}\mathrm{e}}_2{\mathrm{O}}_3+3\mathrm{C}{\mathrm{H}}_4\rightleftarrows$ $6{\mathrm{H}}_2+3\mathrm{C}\mathrm{O}+2\mathrm{F}\mathrm{e}$	252.16	−0.285	298–1300	884.8	[1] *	[59]
4	$2{\mathrm{F}\mathrm{e}}_2{\mathrm{O}}_3+3\mathrm{C}{\mathrm{H}}_4\rightleftarrows$ $6{\mathrm{H}}_2\mathrm{O}+3\mathrm{C}+4\mathrm{F}\mathrm{e}$	158.12	−0.188	298–1300	841.1	[1] *	[59]
5	$3\mathsf{\alpha }-\mathrm{F}\mathrm{e}+\mathrm{C}\rightleftarrows \mathrm{F}{\mathrm{e}}_3\mathrm{C}$	11.98	−0.017	700–1300	704.7	[1] *	analytical solutions for $\mathsf{\alpha }-\mathrm{F}\mathrm{e}$, graphite, ${\mathrm{F}\mathrm{e}}_3\mathrm{C}$ given in [70]
6a 6b	$3\mathsf{\alpha }-\mathrm{F}\mathrm{e}+\mathrm{C}{\mathrm{H}}_4\rightleftarrows {\mathrm{F}\mathrm{e}}_3\mathrm{C}+2{\mathrm{H}}_2$ $3\mathsf{\gamma }-\mathrm{F}\mathrm{e}+\mathrm{C}{\mathrm{H}}_4\rightleftarrows {\mathrm{F}\mathrm{e}}_3\mathrm{C}+2{\mathrm{H}}_2$	101.29 101.29	−0.122 −0.126	773–1273 773–1273	810.4 804.0	[1] [1]	for ${\mathrm{F}\mathrm{e}}_3\mathrm{C}$ assumed after [1] for 485 < T < 1300 K
7	$\mathsf{\alpha }-\mathrm{F}\mathrm{e}\to \mathsf{\gamma }-\mathrm{F}\mathrm{e}$	22.41	−0.066	700–1300	339.5	[85]	-
8 9	$\mathrm{C}\mathrm{O}+{\mathrm{H}}_2\mathrm{O}\rightleftarrows \mathrm{C}{\mathrm{O}}_2+{\mathrm{H}}_2$ ${\mathrm{H}}_2+{\mathrm{C}\mathrm{O}}_2\rightleftarrows \mathrm{C}\mathrm{O}+{\mathrm{H}}_2\mathrm{O}$	−38.06 38.06	0.036 −0.036	298–1300 298–1300	1057.2 1057.2	[1] [1]	Water-Gas-Shift (WGS) [90,91,92] and Reverse-WGS (R-WGS) [90,91,92]
10	$4{\mathrm{H}}_2+\mathrm{C}{\mathrm{O}}_2\rightleftarrows \mathrm{C}{\mathrm{H}}_4+2{\mathrm{H}}_2\mathrm{O}$	−179.56	0.210	298–1300	855.0	[1]	R-WGS as methanation [91,93] Sabatier reaction
11	$\mathrm{C}+{\mathrm{H}}_2\mathrm{O}\rightleftarrows \mathrm{C}\mathrm{O}+{\mathrm{H}}_2$	137.20	−0.156	298–1300, 700–1300 for $\mathrm{C}$	879.5	[1]	WG [4]
12	$\mathrm{C}+\mathrm{C}{\mathrm{O}}_2\rightleftarrows 2\mathrm{C}\mathrm{O}$	172.77	−0.169	298–1400	1022.3	[94]	Bell–Boudouard reaction

* polymorphic phase change No. 7 was assumed.

below summarises stoichiometric reactions and Gibbs free energy equations for reaction (1), as well as multiple correlations for possible polymorphic transformations of iron. For example, positions 2–4 refer to the reaction of ferric oxide with methane according to the equations given in [59], while positions 5–6 are associated with the formation and decomposition of iron carbides (mainly cementite), while position 7 represents the value of the free energy of the $\mathsf{\alpha }-\mathrm{F}\mathrm{e}\to \mathsf{\gamma }-\mathrm{F}\mathrm{e}$ transformation. Analytical approaches to the carburisation and decarburisation reactions of iron have been presented in [62], these results are in agreement with [68] and indicate that the diffusion coefficients of carbon into iron up to 912 °C follow the sequence: $\mathsf{\alpha }-\mathrm{F}\mathrm{e}>\mathsf{\gamma }-\mathrm{F}\mathrm{e}\gg \mathrm{F}{\mathrm{e}}_5{\mathrm{C}}_2$. Interestingly, for the reaction of $\mathrm{F}\mathrm{e}\mathrm{O}$ with methane, in the temperature range of 700–1100/1200 K, the thermodynamic existence of the $\mathrm{F}{\mathrm{e}}_3\mathrm{C}$ phase was determined in [89], while works [66,69,72] discuss the formation of cementite on an active intermediate phase, e.g., on hypereutectoid steels [72] and on orthorhombic crystal structure [60].

The formation of cementite is an endothermic reaction, i.e., its decomposition is accompanied by an exothermic effect. Therefore, for isothermal decomposition of methane, knowledge of the temperature alone seems to be insufficient, while the advantages of following the complete heat balance of the system become more pronounced. Furthermore, when using iron oxide catalysts, the presence of CO₂, CO and water vapour in the reacting system must also be taken into account. For this reason, positions 8–12 in Table 1 present possible pathways involving these species.

A graphical representation of the dependencies obtained from the data presented in Table 1 is given in Figure 2 below. Within the temperature range considered for methane pyrolysis using iron as a catalyst, the homogeneous gas phase reactions such as WGS (No. 8), R-WGS (No. 9) and methanation (No. 10) were omitted. The reason for omitting the WGS reaction (No. 8) is the fact that it can be carried out effectively even at temperatures < 500 °C, while to carry out the methanation reaction (No. 10) at ${∆}_{r}G>0$, a stoichiometric amount of hydrogen is required (specifically ${\mathrm{H}}_2:\mathrm{C}{\mathrm{O}}_2=4:1$), and these conditions do not correspond to the conditions of methane pyrolysis [90,93]. In addition, at temperatures above 950 °C the presence of methane in the reaction atmosphere counteracts methanation. Finally, at $T>T_{eq.}=1091 \mathrm{K} \left(818 °\mathrm{C}\right)$, the ${∆}_{r}G$ of the R-WGS reaction (No. 9) is $<0$, and the equilibrium conversion of $\mathrm{C}{\mathrm{O}}_2$ is $>0.5$ [91]. However, the main reason is that R-WGS is carried out at lower temperatures and requires a bimetallic catalyst [90,91,93]. Similarly to reactions 8–10, reactions 11–12 were omitted from Figure 2 because their occurrence in steady-state conditions is highly unlikely due to the lack of availability of oxygen in the reacting system.

The most likely reaction mechanism (1) suggests that the process takes place at the catalyst surface, specifically by carbon diffusion into the active centre of the catalyst, which for iron catalysts leads to the formation of cementite [4,68]. Looking at the position of the equilibrium temperature ($T_{eq.}$), it can also be observed that iron oxides shift the equilibrium to the right, while metallic iron shifts it slightly to the left. An extreme and very undesirable example of methane decomposition is shown by reaction No. 4 in Table 1, where hydrogen is oxidised to water [59].

3.1. Modelling of the Forward and Backward Reactions

Utilising the phenomenological thermodynamics and in particular model of the first-order kinetics, the equilibrium constant can be expressed in the following logarithmic form:

$$\ln K=\ln k_1-\ln k_{-1}$$

(15)

After multiplying by $\left(-RT\right)$ and using the relation that the free energy of a reaction can be expressed as the difference between the activation energies of the forward and backward reactions [95], we can write that:

$${∆}_r{G}^{\mathrm{\varnothing }}=∆G^{+}-{∆}_{-1}{G}^{+}$$

(16)

Each component of the right-hand side of Equation (16) is obtained by comparing the Arrhenius kinetic constants with Eyring’s TST [96]. Therefore, using the concept presented in [97], we obtain:

$$∆G^{+}=R·T\mathit{ln}{\displaystyle \frac{BT}{k_1}}\mathrm{and} {∆}_{-1}{G}^{+}=R·T\mathit{ln}{\displaystyle \frac{BT}{k_{-1}}}$$

(17)

Substituting Equation (17) into Equation (16) gives the following convenient form of equation:

$${∆}_r{G}^{\varnothing }=R·T\mathit{ln}{\displaystyle \frac{k_{-1}}{k_1}}$$

(18)

These formulae give us a starting point for a qualitative analysis of various processes that show similarities to methane pyrolysis, although it is important to note the following limitations. It must be stressed that this approach is only fully valid in the case of perfectly reversible reactions, i.e., where the reverse reaction yields the same substrate (chemically and morphologically). As already mentioned for the catalytic process of methane decomposition, it can only be applied here with due caution, since the reactions taking place on the surface of the catalyst determine its course (e.g., through the structure of the carbon deposits formed and their reactions with the catalyst itself) and it is virtually impossible to fulfil the reversibility conditions completely. Further insight into this issue can be found, for example, in [98], where the importance of the third term on the right-hand side of Equation (16) ${∆}_n{G}^{\varnothing }$ being associated with solid phase transformations, such as nucleation, is adequately analysed. For the sake of the following discussion, here the term ${∆}_n{G}^{\varnothing }$ was tentatively assumed to be zero. However, as will be shown in the following paragraphs, the error introduced by such an approximation is not remarkable.

3.2. Variant 1—From Iron Carburization to Methane Decomposition on the Iron Surface

According to [62], the above-mentioned processes of carburisation and decarburisation of iron are based on isothermal reactions carried out in the temperature range of 800–1040 °C. In the original text, the authors use the following notation:

$$\mathrm{C} \left(\mathrm{i}\mathrm{n} \mathsf{\gamma }-\mathrm{F}\mathrm{e}\right) +2{\mathrm{H}}_2\left(\mathrm{g}\right) \leftrightarrow \mathrm{C}{\mathrm{H}}_4\left(\mathrm{g}\right)$$

(19)

The enthalpy of this reaction was determined to be ${∆}_r{H}^{\varnothing }=-124.3$ kJ·mol⁻¹. In [62] the rate constants for this reaction were also determined, but it is worth noting that when cross-referencing with this review, attention must be paid to the differences in notation and units. In the following considerations, the reverse direction of reaction (19) is used, i.e., from right to left. Furthermore, by unifying the dimensions to the scale ${\mathrm{s}}^{-1}$, we obtain the following rate constant for the reaction of methane decomposition:

$$k_1=4.55·{10}^3 \left[\exp \left(-{\displaystyle \frac{\mathrm{230,120}}{\mathrm{R}T}}\right)\right],$$

while for the reaction of methane synthesis [as indicated by the direction of Equation (19)]:

$$k_{-1}=2.6·{10}^{-2}\left[\exp \left(-{\displaystyle \frac{\mathrm{105,860}}{\mathrm{R}T}}\right)\right]$$

Substitution of these two reaction rate constants into Equation (18) yields:

$${∆}_r{G}^{\mathrm{\varnothing }}=124.27-0.100T$$

(20)

Interestingly, the parameters of Equation (20) obtained in this way are similar to the values given for position 6b in Table 1. For ${\mathrm{F}\mathrm{e}}_3\mathrm{C}$ the slightly lower value of the intercept can be explained by the non-linear relationship between free energy and temperature [70], which in a simplified approach for the range of 298.15–1500 K can be approximated by a second-degree parabolic function [1,60]. For smaller temperature ranges the polynomial can be reduced to a linear form and thus the values of the intercept can be determined to vary in the range 70–116 $\mathrm{k}\mathrm{J}·{\mathrm{m}\mathrm{o}\mathrm{l}}^{-1}$. If, on the other hand, the 500–2000 K range is accepted, the free energy equation for position 6b given in Table 1 takes the form: ${∆}_r{G}^{\varnothing }=126.27-0.144·T$. Therefore, it must be emphasised that the above expression of Equation (20) correlates well with the reaction of methane decomposition in the presence of metallic iron catalyst ($\mathsf{\gamma }-\mathrm{F}\mathrm{e}$). Unfortunately, the authors in [60] do not specify to what extent it is possible to treat the solution of $\mathsf{\gamma }-\mathrm{F}\mathrm{e}$ and $\mathrm{C}$ as cementite.

Finally, returning to the notation of Equation (19), the above considerations show that even with the correction for the reversibility of reaction (1), the correctness of analysing the process of methane decomposition using Eyring’s TST [96] still holds. In this context, the above statement is consistent with the basic principles of reaction (1) given by [7].

3.3. Variant 2—Free Energy of Activation of Reaction (1)

The above analysis shows the link between the work carried out to develop the carburisation/decarburisation process and the methane pyrolysis process. For this reason, it is now important to contrast this similarity with the results of an isolated analysis of the general scheme of methane pyrolysis, as derived from the reaction notation (1).

Similarly to the previous case, in the following example two different sources of kinetic data, free of the potential influence of Fe as catalyst, are analysed for their convergence. Because formally, now the value of the component $∆G^{+}$ is sought, thus the notification of Equation (16) changes to:

$$∆G^{+}={∆}_r{G}^{\mathrm{\varnothing }}+{∆}_{-1}{G}^{+}$$

(21)

As the first example, data from [99] can be used. Here the authors described a process of high temperature (1080–3000 °C) and high pressure (up to 5.61 MPa) methane synthesis. When hydrogen is reacted with graphitic materials the following kinetic parameters were obtained $A=5.9 {\mathrm{s}}^{-1}$, $E=101.50 \mathrm{k}\mathrm{J}·\mathrm{m}\mathrm{o}{\mathrm{l}}^{-1}$. If we restrict the temperature range close to the region of our interests, i.e., 500–1500 K, these data give the following expression for the free energy of activation:

$${∆}_{-1}{G}^{+}=93.88+0.248·T$$

(22)

Substitution of the above expression and the data from Table 1, position 1, into Equation (21) gives:

$$∆G^{+}=181.92+0.140·T$$

(23)

A second example can be found in [55]. In this work, the authors have determined and used data describing the kinetics of methane decomposition using a graphite crucible and a liquid slag with a predominant content of CaO (43% by mass) and only a minimal content of FeO (0.37% by mass). This work provides the following kinetic parameters, $A=4.053·{10}^6 \mathrm{c}{\mathrm{m}}^{-1}·{\mathrm{s}}^{-1}$, $E=190.31 \mathrm{k}\mathrm{J}·\mathrm{m}\mathrm{o}{\mathrm{l}}^{-1}$, which can thus be used to determine $∆G^{+}$:

$$∆G^{+}=182.69+0.130·T$$

(24a)

To keep the notation of the units consistent, the data must be corrected to the form $A=4.819·{10}^7 {\mathrm{s}}^{-1}$ and as a result, we obtain:

$$∆G^{+}=181.41+0.120·T$$

(24b)

Although the input data used to carry out variant 2 of our analysis came from two very different experimental approaches, they both focused on the same reaction principle (1). Both gave not only similar coefficients, but also coefficients close to those previously determined for the carburisation of iron. Moreover, despite the fact that the references [62,99] are related to iron metallurgy and therefore not directly aligned with the main scope of this review, based on the kinetic data provided, it is possible to consider the solid of $\mathrm{F}\mathrm{e}+\mathrm{C}$ as a catalyst for methane decomposition.

The demonstrated equivalence of analyses of data derived from different sources, together with the results of their phenomenological thermodynamic assessment, represents a significant step towards using the knowledge and insights gained so far in the context of kinetic considerations.

4. Kinetics—Introduction

In the context of kinetics, the literature on methane pyrolysis includes such variables as: the methane and hydrogen mole fractions, their partial pressures, the methane conversion rate, as well as the rate (or amount) of carbon formation per unit of catalyst mass—expressed in [$\mathrm{g}$ or $\mathrm{m}\mathrm{o}\mathrm{l}·{\mathrm{g}}_{\mathrm{c}\mathrm{a}\mathrm{t}}^{-1}·{\mathrm{u}\mathrm{n}\mathrm{i}\mathrm{t} \mathrm{o}\mathrm{f} \mathrm{t}\mathrm{i}\mathrm{m}\mathrm{e}}^{-1}$] or [$\mathrm{g}$ or $\mathrm{m}\mathrm{o}\mathrm{l}·{\mathrm{g}}_{\mathrm{c}\mathrm{a}\mathrm{t}}^{-1}$], respectively. The relationship between the amount of carbon deposit formed is derived from a volumetric balance, where for reaction (1) $\dot{V}= const$, so the unit typically used in the literature is [$\mathrm{d}{\mathrm{m}}^3$ or ${\mathrm{m}}^3$·${\mathrm{g}}_{\mathrm{c}\mathrm{a}\mathrm{t}}^{-1}·{\mathrm{u}\mathrm{n}\mathrm{i}\mathrm{t} \mathrm{o}\mathrm{f} \mathrm{t}\mathrm{i}\mathrm{m}\mathrm{e}}^{-1}$]. For further considerations, it should be noted that parameters such as the degree of conversion should be considered particularly important due to their practical nature. Following the analyses presented in [27,58], we can integrate the approach to quantify the carbon deposit by using the following term:

$$\left[C\right]={\displaystyle \frac{12}{22.4}}\dot{V}{\int }_0^t\alpha \mathrm{d}t$$

(25)

As a different approach to modelling the kinetics of carbon deposition, the following expression can be found in the literature [35,39,42,100,101,102]:

$$r_C={\displaystyle \frac{\mathrm{d}\left[C\right]}{\mathrm{d}t}}=k(P_{\mathrm{C}{\mathrm{H}}_4}{)}^n,n=0.5÷1,$$

(26)

Here $r_C$ can be understood as the initial rate of formation of the carbon deposit. This approach is very useful because it depends on only one variable ($P_{\mathrm{C}{\mathrm{H}}_4}$). Considering the research presented in [101], it follows that ${\displaystyle \frac{P_{\mathrm{C}{\mathrm{H}}_4}}{P^{\mathrm{\varnothing }}}}=x_0=\mathrm{c}\mathrm{o}\mathrm{n}\mathrm{s}\mathrm{t}$ for a single methane carrier, but when we are dealing with several different methane streams $x_0=\mathrm{v}\mathrm{a}\mathrm{r}$, and this complicates the analysis.

Furthermore, in works [14,83], Equation (26) was used to indicate the existence of a maximum rate of formation of the carbon deposit. The value of this parameter is a direct consequence of the activity of the catalyst used. Once the maximum is reached, a decrease in the rate of formation of the $C$ deposit is observed; therefore, as indicated in [103], the formation of carbon can also be expressed empirically as a function of time:

$$\left[C\right]\propto t^n, \mathrm{n}<1$$

(27)

A better understanding of this relationship can be gained by substituting the simplest linear relationship for the conversion degree versus time in Equation (25). When integrated, the relationship between the amount of $C$ deposited over time takes the form of a second-order parabolic equation. Importantly, however, from a practical point of view, another useful way of presenting Equation (25) is to replace the integral with the average conversion degree. In this case, we can work with the following simple expression:

$$\left[C\right]={\displaystyle \frac{12}{22.4}}\dot{V}\overline{\alpha }t$$

(28)

Table 2. Comparison of experimental data for fluidised bed pyrolysis of methane using NiCu/Al₂O₃ as catalyst.

Temperature in °C; Model for H2	$\overline{\mathit{\alpha }}$, Calculated as an Average Value of an Integral	$\left[\mathit{C}\right], \mathbf{g}·{\mathbf{g}}_{\mathbf{c}\mathbf{a}\mathbf{t}}^{-1}$ Calculated from Equation (28)	$\left[\mathit{C}\right],·{\mathbf{g}·\mathbf{g}}_{\mathbf{c}\mathbf{a}\mathbf{t}}^{-1}$ in [27] (Time of the Process in Minutes)
650; exponential	0.127	5.71	4.97 (420)
700; linear	0.404	18.17	18.41 (420)
750; exponential	0.288	12.95	12.95 (420)
800; exponential	0.303	8.44	8.80 (300)

below summarises a selection of literature data on the amount of $C$ deposited during the fluidised bed methane pyrolysis process [27]. In these tests, the process was run at temperatures of 650–800 °C using NiCu/Al₂O₃ catalyst and with $WHSV=12 \mathrm{L}·{\mathrm{g}}_{\mathrm{c}\mathrm{a}\mathrm{t}}^{-1}·{\mathrm{h}}^{-1}$.

The second column of Table 2 shows the calculated integral mean for methane conversion, integrated from time $t=0$, when the point of maximum conversion is reached, to the end of the experiment (finite time). The third and fourth columns show, respectively, the amounts of carbon deposition calculated using Equation (28) and the experimental amounts of deposition reported in [27]. These data refer directly to the phase of methane pyrolysis when the catalyst is deactivated (the second phase of the process) and therefore do not take into account the initial moments of reaction (1).

In contrast to the example given above the literature provides descriptions of methane pyrolysis processes where both the initial and deactivation steps are considered. During the initial stage, methane conversion, carbon deposition and hydrogen content in the process gas increase while the methane concentration decreases. This period can therefore be considered as the activation step. While during the second period of the process, after time $t\ge {\tau }_D$, where ${\tau }_D$ is the deactivation time, a continuous process of catalyst deactivation is observed [104]. The time dependence of the reaction rate for the deactivation period can be described, for example, as follows

$$r\left(t\right)=r_0\mathrm{e}\mathrm{x}\mathrm{p}\left(-{\displaystyle \frac{t}{{\tau }_D}}\right), \mathrm{for} t=0,r\left(t\right)=r_0$$

(29)

Following after [105] or [106], which state that for many process conditions, the influence of the reverse reaction is negligible and can therefore be omitted, we can also derive the following expression for the methane conversion rate:

$$r=-{\displaystyle \frac{\mathrm{d}{x}_{\mathrm{C}{\mathrm{H}}_4}}{\mathrm{d}t}}=k_1(x_{\mathrm{C}{\mathrm{H}}_4}{)}^n, n=1÷1.5$$

(30)

The inadequacy of the above attempts to describe the process of methane pyrolysis calls for an extension of the discourse to include the concept of elementary reaction.

4.1. Kinetics of the Elementary Reaction (1)

The IUPAC Gold Book defines an elementary reaction as follows [107]: “A reaction for which no reaction intermediates were detected or need to be postulated in order to describe the chemical reaction on a molecular scale. An elementary reaction is assumed to occur in a single step and to pass through a single transition state”.

In practice, the above description can be restricted to the complete conversion of substrates according to the stoichiometric expression or to a limited equilibrium state.

In order to determine whether reaction (1) can be treated as an elementary reaction and what its kinetic orders are when considering both limits of the methane conversion degree, the following dimensionless notations were used in the following discussion: $x_{\mathrm{C}{\mathrm{H}}_4}$—molar concentration of methane,

$x_{\mathrm{C}{\mathrm{H}}_4\_0}$—molar concentration of methane at $t=0$,

$x_{\mathrm{C}{\mathrm{H}}_4\_eq}$—equilibrium molar concentration of methane,

$\alpha$—conversion degree of methane,

${\alpha }_{eq}$—equilibrium conversion degree of methane,

as well as,

$k_1$—forwards reaction rate constant ${\mathrm{s}}^{-1}$,

$k_{-1}$—backwards reaction rate constant ${\mathrm{s}}^{-1}$,

$k_3$—reaction rate constant of a different kind, ${\mathrm{s}}^{-1}$.

For reaction (1) the equilibrium constant is determined from the relationship characteristic of first-order kinetics (see Equation (15)):

$$K={\displaystyle \frac{k_1}{k_{-1}}}$$

(31)

The rate of such a reaction can be expressed in the following form:

$$r=k_1{x}_{\mathrm{C}{\mathrm{H}}_4}-k_{-1}(1-x_{\mathrm{C}{\mathrm{H}}_4}{)}^2, T=\mathrm{c}\mathrm{o}\mathrm{n}\mathrm{s}\mathrm{t}$$

(32)

where:

$$r=-{\displaystyle \frac{\mathrm{d}{x}_{\mathrm{C}{\mathrm{H}}_4}}{\mathrm{d}t}}$$

(33)

The coupling of Equations (31)–(33) leads to the form of an equation with separate variables:

$$-{\displaystyle \frac{\mathrm{d}{x}_{\mathrm{C}{\mathrm{H}}_4}}{\mathrm{d}t}}=k_1\left[x_{\mathrm{C}{\mathrm{H}}_4}-{\displaystyle \frac{(1-x_{\mathrm{C}{\mathrm{H}}_4}{)}^2}{K}}\right]$$

(34)

In Equation (34) the equilibrium constant $K$ can be expressed eg. using Equation (6) or (8).

Importantly for further considerations in the work of [82] using the previous publication [108], the authors have carried our research in a manner most similar to the analysis of an elementary reaction, thus a semi-batch system was used, the methane feed was kept minimal, only to maintain constant pressure, and the hydrogen concentration measurements were taken after the process. In this way, for the temperature from 700 °C to 900 °C and total pressure of 2.82–5.61 MPa, these results can be attributed to a gradient-free process of methane pyrolysis. In this work, the following initial condition was accepted $x_{\mathrm{C}{\mathrm{H}}_4\_0}$:

$${g\left(x_{\mathrm{C}{\mathrm{H}}_4}\right)=k}_{1}t$$

(35)

and thus for $t=0, x_{\mathrm{C}{\mathrm{H}}_4}=x_{\mathrm{C}{\mathrm{H}}_4\_0} \mathrm{a}\mathrm{n}\mathrm{d} x_{\mathrm{C}{\mathrm{H}}_4}<x_{\mathrm{C}{\mathrm{H}}_4\_eq}$ (intercept = 0), the following expression was determined:

$$k_{1}t={\displaystyle \frac{x_{\mathrm{C}{\mathrm{H}}_4\_0}-x_{\mathrm{C}{\mathrm{H}}_4\_eq}}{x_{\mathrm{C}{\mathrm{H}}_4\_0}+x_{\mathrm{C}{\mathrm{H}}_4\_eq}}}\ln \left[{\displaystyle \frac{{x_{\mathrm{C}{\mathrm{H}}_4\_0}}^2-x_{\mathrm{C}{\mathrm{H}}_4\_eq}}{x_{\mathrm{C}{\mathrm{H}}_4\_0}\left(x_{\mathrm{C}{\mathrm{H}}_4}-x_{\mathrm{C}{\mathrm{H}}_4\_eq}\right)}}\right]$$

(36)

Now by substituting the relationship between the molar concentration of methane and its degree of conversion:

$$x=x_{\mathrm{C}{\mathrm{H}}_4\_0}{\displaystyle \frac{1-\alpha }{1+x_{\mathrm{C}{\mathrm{H}}_4\_0}\alpha }},$$

(37)

into Equation (36), we obtain the degree of methane conversion that is related to the kinetic form and thus the following can be obtained:

$$k_{1}t={\displaystyle \frac{{\alpha }_{eq}\left(1+x_{\mathrm{C}{\mathrm{H}}_4\_0}\right)}{2-{\alpha }_{eq}\left(1-x_{\mathrm{C}{\mathrm{H}}_4\_0}\right)}}·\ln \left[{\displaystyle \frac{{\alpha }_{eq}+\alpha -\alpha {·\alpha }_{eq}\left(1-x_{\mathrm{C}{\mathrm{H}}_4\_0}\right)}{{\alpha }_{eq}-\alpha }}\right]$$

(38)

For $x_{\mathrm{C}{\mathrm{H}}_4\_0}=1$ (pure methane) Equation (38) becomes:

$$k_{1}t={\alpha }_{eq}\ln \left({\displaystyle \frac{{\alpha }_{eq}+\alpha }{{\alpha }_{eq}-\alpha }}\right),{\alpha }_{eq}>\alpha ,$$

(39)

and, given that ${\displaystyle \frac{\alpha }{{\alpha }_{eq}}}=z$, we can alternatively write Equation (39) as:

$${k_{1}t=\alpha }_{eq}\ln \left({\displaystyle \frac{1+z}{1-z}}\right),z<1$$

(40)

Knowing that $z$ corresponds to the thermodynamic efficiency of reaction (1), in the following discussion we can use the kinetic functions $g\left(\alpha \right)$ and $g\left(z\right)$ to relate to the notation of their complete kinetic forms ($=k_{1}t$).

In Equation (40) for $t=0, z=0$, or for $t\to \infty , z\to 1,$ and thus for small values of $z$, an alternative solution can be sought by expanding to a series:

$$\ln \left({\displaystyle \frac{1+z}{1-z}}\right)=2\sum _1^{\mathrm{\infty }}\left[\left({\displaystyle \frac{1}{2n-1}}\right)z^{2n-1}\right],z^2<1$$

(41)

Furthermore, an approximation to Equation (40) can be found by limiting this series to its first term:

$$\ln \left({\displaystyle \frac{1+z}{1-z}}\right)\cong 2z$$

(42)

and by substituting Equation (42) into Equation (40), the following expression can be used as a way of describing a minimal operating line:

$$\alpha ={\displaystyle \frac{k_1}{2}}t=k_{z}t$$

(43)

The following Figure 3 can be used as an example to describe the effectiveness of this procedure. The figure shows a transformation of the experimental results presented in [82,108] to the form given by Equation (39) and its further simplification to the linear relation of Equation The result of this exercise is that for the conditions: $x_{\mathrm{C}{\mathrm{H}}_4\_0}=1$, methane flow $\dot{V}\approx 0$, then $\alpha =1-x_{\mathrm{C}{\mathrm{H}}_4},$ the kinetic parameters thus calculated are: $E=130.5$ kJ·mol⁻¹; $\ln A=8.60$ ($A$ in ${\mathrm{s}}^{-1}$). For reference, according to [82,108], the kinetic parameters of methane pyrolysis are: $E=131.0$ kJ·mol⁻¹; and $\ln A=8.59$ ($A$ in ${\mathrm{s}}^{-1}$). This can be used as evidence that in the context of the methane pyrolysis process, the above-described simplified approach to kinetic modelling can produce very accurate results.

The commentary on the data presented in Figure 3 can be further summarised in a few remarks. Comparing the structure of Equations (36) and (43), it is clear that the kinetics of the process can be interpreted using several models of varying complexity. This is a fairly common observation in studies of heterogeneous systems involving both solid and gas phases. For example, if Equation (34) is expressed in another form, say as:

$$-{\displaystyle \frac{\mathrm{d}{x}_{\mathrm{C}{\mathrm{H}}_4}}{\mathrm{d}t}}=k_1{x}_{\mathrm{C}{\mathrm{H}}_4}\left[1-{\displaystyle \frac{(1-x_{\mathrm{C}{\mathrm{H}}_4}{)}^2}{x_{\mathrm{C}{\mathrm{H}}_4}K}}\right]$$

(44)

Using Equations (4) and (9), we can rewrite Equation (44) as:

$$-{\displaystyle \frac{\mathrm{d}{x}_{\mathrm{C}{\mathrm{H}}_4}}{\mathrm{d}t}}=k_1{x}_{\mathrm{C}{\mathrm{H}}_4}\left(1-{\displaystyle \frac{K_x}{K}}\right)$$

(45)

If reaction (1) remains “very far” from equilibrium $\left(K\gg K_x\right)$ according to [107], it follows first-order kinetics for $\dot{V}=0$, and thus in terms of the degree of methane conversion, Equation (45) can be described as:

$$g\left(\alpha \right)=-\ln \left(1-\alpha \right)=k_{1}t,T=const$$

(46)

On the other hand, this also shows that reaction (1) has zero order kinetics when the degree of conversion of methane is very small.

In kinetic analyses, it is usually considered appropriate to accept models with the least complexity while maintaining sufficient accuracy. The results of studies of the elementary reaction (1) under high pressure [82,108] indicate that the equilibrium state, represented by the equilibrium conversion degree, is the limiting factor here. From the fact that the deviation from equilibrium influences the kinetics, we know that such processes can be expressed by different modelling approaches, e.g., according to the $n^{th}$ or fractional order variability [105].

As mentioned above, the literature on the topic of methane pyrolysis refers to the inequality $z={\displaystyle \frac{\alpha }{{\alpha }_{eq.}}}<1$ to describe the thermodynamic efficiency of the process, but for certain process conditions, e.g., for some temperature ranges, $z\approx 1$ [19,38]. For the research published in [38] such a case was found, e.g., for the catalysts consisting of 50% Ni-10% Cu/SiO₂.

4.2. Continuous Process—Conditions of Constant Flow of Methane

For the sake of broadening the discussion to continuous flow systems, it is also needed to find a solution to Equation (34) when the flow of gases through the reactor are non-zero, $\dot{V}>0$. We can achieve this by differentiating Equation (37) with respect to time, which reduces it to the expression for the rate of reaction (1):

$$-{\displaystyle \frac{\mathrm{d}{x}_{\mathrm{C}{\mathrm{H}}_4}}{\mathrm{d}t}}=x_{\mathrm{C}{\mathrm{H}}_4\_0}\left(1+x_{\mathrm{C}{\mathrm{H}}_4\_0}\right)·(1+x_{\mathrm{C}{\mathrm{H}}_4\_0}\alpha {)}^{-2}{\displaystyle \frac{\mathrm{d}\alpha }{\mathrm{d}t}}$$

(47)

After substituting Equations (37) and (47) as well as the thermodynamic equilibrium constant given by Equation (12) into Equation (34), when the condition $z={\displaystyle \frac{\alpha }{{\alpha }_{eq}}}<1$ is fulfilled, we finally obtain the following kinetic expression for reaction (1) under flow conditions in terms of the degree of methane conversion:

$${\displaystyle \frac{\mathrm{d}\alpha }{\mathrm{d}t}}=k_1{\displaystyle \frac{2-{\alpha }_{eq}\left(1-x_{\mathrm{C}{\mathrm{H}}_4\_0}\right)}{{\alpha }_{eq}\left(1+x_{\mathrm{C}{\mathrm{H}}_4\_0}\right)}}·{\displaystyle \frac{\left(c·\alpha +{\alpha }_{eq}\right)\left({\alpha }_{eq}-\alpha \right)}{{\alpha }_{eq}\left(1+c\right)}}, \mathrm{where} c=1-{\alpha }_{eq}\left(1-x_{\mathrm{C}{\mathrm{H}}_4\_0}\right)$$

(48)

For the initial condition $x_{\mathrm{C}{\mathrm{H}}_4\_0}=1$ (pure methane), the complexity of Equation (48) is reduced to the following, much simpler form for interpretation:

$${\displaystyle \frac{\mathrm{d}\alpha }{\mathrm{d}t}}={\displaystyle \frac{{\mathrm{k}}_1}{2}}\left[1-{\left({\displaystyle \frac{\mathsf{\alpha }}{{\mathsf{\alpha }}_{\mathrm{e}\mathrm{q}}}}\right)}^2\right], T=\mathrm{c}\mathrm{o}\mathrm{n}\mathrm{s}\mathrm{t}$$

(49)

which can be similarly expressed in terms of thermodynamic efficiency:

$${\displaystyle \frac{\mathrm{d}z}{\mathrm{d}t}}={\displaystyle \frac{k_1}{2{\alpha }_{eq}}}\left(1-z^2\right)$$

(50)

Equations (49) and (50) indicate that under flow conditions of methane, the reaction rate is quasi-second order, which is consistent with the exponent related to the hydrogen formation presented in Equation (32). Noteworthy, the considerations expressed through Equations (47)–(50) highlight that depending on the conditions under which reaction (1) is carried out, its order can vary.

The solution to equation (50) can be sought in a similar way to that given for Equation (40). The form of Equation (40) raises the question of whether this expression is suitable for use in complex studies of reaction (1), especially where the presence of catalysts and their interactions with methane are taken into account.

On the other hand, for $({\alpha }_{\mathit{eq}}=1)$ and $(x_{\mathrm{C}{\mathrm{H}}_4\_0}<1)$, after substitution and rearrangement of Equation (38), we obtain:

$$g\left(\alpha \right)=k_{1}t=\ln \left({\displaystyle \frac{1+x_{\mathrm{C}{\mathrm{H}}_4\_0} \alpha }{1-\alpha }}\right)$$

(51)

Therefore, for $t=0, \alpha =0$, and for $t\to \mathrm{\infty }, \alpha \to 1$, which is consistent with the results obtained for the non-flow conditions and the thermodynamic efficiency of reaction (1). This leads us to the possibility of further extending the range of conditions for which we can analyse the course of reaction (1), and so in the following section the issue of using methane diluted with inert carrier gas as a feed will be analysed

4.3. Continuous Process—Dilute-Methane Feed Conditions

The reasoning for such an experimental approach derived directly from stoichiometry and Le Chatelier’s principle. This problem can be tackled using a similar analytical approach to the one presented above for the process flow conditions. Here taking into account the mass balance of Equation (37) and using Equation (12), the following expression for the equilibrium constant related to the initial molar fraction of methane $x_{\mathrm{C}{\mathrm{H}}_4\_0}<1$ can be derived:

$$K={\displaystyle \frac{(1-x_{\mathrm{C}{\mathrm{H}}_4\_0}+2x_{\mathrm{C}{\mathrm{H}}_4\_0}{\alpha }_{eq}{)}^2}{x_{\mathrm{C}{\mathrm{H}}_4\_0}\left(1+x_{\mathrm{C}{\mathrm{H}}_4\_0}{\alpha }_{eq}\right)\left(1-{\alpha }_{eq}\right)}}$$

(52)

From Equation (52) equilibrium conversion degree reduces to the following quadratic equation:

$${\alpha }_{eq}^2+{\displaystyle \frac{1-x_{\mathrm{C}{\mathrm{H}}_4\_0}}{x_{\mathrm{C}{\mathrm{H}}_4\_0}}}{\alpha }_{eq}+{\displaystyle \frac{{x_{\mathrm{C}{\mathrm{H}}_4\_0}}^2-x_{\mathrm{C}{\mathrm{H}}_4\_0}\left(2+K\right)+1}{{x_{\mathrm{C}{\mathrm{H}}_4\_0}}^2\left(4+K\right)}}=0$$

(53)

As an example, using Equation (53) it can be determined that for 50% dilution of methane with nitrogen $x_{\mathrm{C}{\mathrm{H}}_4\_0}=0.5$ and for temperature of $T=1173 \mathrm{K}$, we can calculate from Equation (8) that $K=K_p=55.40$. Using this value for the equilibrium constant, the following expression can be given:

$${\alpha }_{eq}^2+{\alpha }_{eq}-1.848=0$$

(54)

After solving we obtain that ${\alpha }_{eq}=0.950\cong 1$. This result showcases that it is justified to analyse the kinetics of methane pyrolysis while neglecting the effect of the reversibility of reaction (1), i.e., setting the term $k_{-1}=0$. Besides the dilution of the methane stream, also an increase in temperature $\left(K\to \mathrm{\infty }\right)$ justifies such an analytical approach.

Accepting the condition of the irreversibility of reaction (1), by combining Equations (32), (33), (37) and (47) the following kinetic equation can be obtained:

$${\displaystyle \frac{\mathrm{d}\alpha }{\mathrm{d}t}}={\displaystyle \frac{k_1}{1+x_{\mathrm{C}{\mathrm{H}}_4\_0}}}\left(1-\alpha \right)\left(1+x_{\mathrm{C}{\mathrm{H}}_4\_0}\alpha \right), T=\mathrm{c}\mathrm{o}\mathrm{n}\mathrm{s}\mathrm{t}$$

(55)

Equation (55) is a solution to Equation (51) and for $x_{\mathrm{C}{\mathrm{H}}_4\_0}=1$ and ${\alpha }_{eq}=1$ solutions for Equation (51) are identical.

Analogous consideration for Equation (43) leads to the following expression:

$$\alpha ={\displaystyle \frac{k_1}{1+x_{\mathrm{C}{\mathrm{H}}_4\_0}}}t=k_{z}t$$

(56)

Table 3. Collation of the series of equations derived based on Equation (34).

No	Equation	$\mathbf{Applies} \mathbf{to} \mathbf{Kinetic} \mathbf{Functions} \mathit{g}\left(\mathit{\alpha }\right)$$\mathbf{or} \frac{\mathit{d}\mathit{\alpha }}{\mathit{d}\mathit{t}}$
1	(38)	$\dot{V}=0$$; x_{\mathrm{C}{\mathrm{H}}_4\_0}<1$$; \alpha <{\alpha }_{eq}<1$
2	(39)	Equation (38) for $\dot{V}=0$$; x_{\mathrm{C}{\mathrm{H}}_4\_0}=1$$; \alpha <{\alpha }_{eq}<1$,
3	(51)	Equation (38) for $\dot{V}=0$$; x_{\mathrm{C}{\mathrm{H}}_4\_0}<1;$$\alpha <{\alpha }_{eq}<1$
4	(48)	$\dot{V}>0$$; x_{\mathrm{C}{\mathrm{H}}_4\_0}<1$$; \alpha <{\alpha }_{eq}<1$
5	(49) and (50) after differentiating vs. Equation (48)	$\dot{V}>0$$; x_{\mathrm{C}{\mathrm{H}}_4\_0}=1$$; \alpha <{\alpha }_{eq}<1$
6	(40)	Solution of Equation (50)
7	(55)	$\mathrm{Solution} \mathrm{of} \mathrm{Equation} \left(51\right) \mathrm{for} k_{-1}=0$
8	(41)–(47), (52)–(54), (56)	Equations supporting the discussion

below presents the series of equations derived from the fundamental equation (Equation (34)) and its solution (Equation (36)) in their integrated ($g\left(\alpha \right))$ and differential forms (${\displaystyle \frac{d\alpha }{dt}}$), derived for specific process conditions: flow rate ($\dot{V}$), initial mole fraction of methane ($x_{\mathrm{C}{\mathrm{H}}_4\_0}$), and thermodynamic constraint (${\alpha }_{eq})$.

Table 3 shows the interweaving of different forms of kinetic equations and highlights their common elements depending on the basic process conditions.

The common element sought concerns first-order kinetics in integral form, and the expression $[-\ln \left(1-\alpha \right)]$ does not appear explicitly in Equations (36), (38), (39) and (51). The following

Table 4. Integral forms of $g\left(\alpha \right)=k_{1}t$ given in the literature.

No.	$\mathbf{d}\mathit{\alpha }/\mathbf{d}\mathit{t}$=	$\mathit{g}\left(\mathit{\alpha }\right)=$	Remarks	Sources
1	${\displaystyle \frac{k_1}{1+x_{\mathrm{C}{{\mathrm{H}}_4}_0}}}\left(1-\alpha \right)$ $\left(1+x_{\mathrm{C}{\mathrm{H}}_4\_0}\alpha \right)$	$\ln \left({\displaystyle \frac{1+x_{\mathrm{C}{\mathrm{H}}_4\_0}\alpha }{1-\alpha }}\right)$	For $x_{\mathrm{C}{\mathrm{H}}_4\_0}=1$; $\mathrm{d}\alpha /\mathrm{d}t={\displaystyle \frac{k_1}{2}}\left(1-{\alpha }^2\right)$	Equations (51) and (55)
2	$k_1{\displaystyle \frac{1-\alpha }{1+\alpha }}$	$-2\ln \left(1-\alpha \right)-\alpha$	Proposed for analysis of reactivity of cokes towards CO2; year 1955 [109]	Eq. proposed by Lee [102], model D–J
3	$k_1{\displaystyle \frac{1-\alpha }{1+\epsilon \alpha }}$	$-\left(1+\epsilon \right) \ln \left(1-\alpha \right)-\epsilon \alpha$, $\epsilon =1/2$	1. $\epsilon$ in [102] as well as in [22,110] is defined in the form given by [111]. 2. for system consisting of isomorphic spheres $\epsilon =1-\left(\pi /6\right)=0.476$	[22] (2004) and [110] (2016)
4	$\left({\displaystyle \frac{k_1}{\beta }}\right){\displaystyle \frac{1-\alpha }{1+x_{\mathrm{C}{\mathrm{H}}_4\_0}\alpha }}$	$-\left(1+x_{\mathrm{C}{{\mathrm{H}}_4}_0}\right)$ $\beta \ln \left(1-\alpha \right)-x_{\mathrm{C}{\mathrm{H}}_4\_0}\beta \alpha$	Acc. to [112], $\beta$ = 0.43–0.71	[112,113,114], for $\beta =1$, $T=1023–1373 \mathrm{K}$, for three types of reactors [115,116]
5	$k_1\left(1-\alpha \right)$	$-\ln \left(1-\alpha \right)$	Here $\epsilon =0$ which acc. to [111] means constant density	Equation (46), model F1

presents similar considerations found in the literature. Importantly, these equations refer to the following conditions: $\dot{V}>0$, $x_{\mathrm{C}{\mathrm{H}}_4\_0}<1$, $\alpha \in \left[0;{\alpha }_{max}\right]$, ${\alpha }_{eq}=1$, while positions no. 2–3 consider the catalytic process. For practical reasons, it is important that these expressions include the term $[-\ln \left(1-\alpha \right)]$, which is characteristic of first-order kinetics. In the following considerations, this term will be denoted by the symbol (F1).

The reaction rate shown in Table 4 indicates that the reaction order ranges from first-order kinetics (model F1, No. 5) to quasi-second order (Nos. 1 and 2). The process described in position 4 is particularly interesting as it uses a solar heating system (solar reactor) and is carried out without the use of a catalyst. The uniqueness of this process lies in the injection of methane into argon at high temperatures. The residence time in the reactor was calculated as a surrogate time, which is the ratio of the reactor volume to the flow rate. The physical dilatation factor $\left(\beta \right)$ [112,113,114] takes into account the relationship between the state parameters (temperature, pressure) at the reactor inlet and those under reaction conditions. Further explanations and implications of this approach are analysed in [115,116]. It is important to note that in these works the calculated Arrhenius law parameters are similar in value to previous studies. In addition, the analysis of the equations given in [112,115,116] indicates that the chemical expansion factor directly corresponds to the value described here as ($x_{\mathrm{C}{\mathrm{H}}_4\_0}$). The results of the studies presented in [112,113,114,115,116] are therefore significant enough to be considered of primary importance. Finally, as a result of the considerations presented in Table 4, it can be noted that with ε = 0, position No. 5 and Equation (46) stand as the simplest kinetic model of reaction (1), model (F1). This model stands out from the rest as it characterises a separate kinetic approach and needs to be fitted with additional process-related terms to adequately characterise reaction (1).

The comprehensive nature of the above considerations on the kinetic reaction model (1) provides a basis on which experimental data can be analysed. The following examples illustrate their validity in practice.

The result of the first example is shown in Figure 4. The data for the analysis were extracted from the experimental research presented in [57]. The author proposed the use of a specially prepared iron mesh for the catalytic decomposition of methane. The pyrolysis was performed at $t=800 °\mathrm{C}, { x}_{\mathrm{C}{\mathrm{H}}_4\_0}=1$; however, the catalyst was previously activated with hydrogen. Analysis of the data set using the approach presented in Figure 3, which is only relevant to the elementary reaction, gives unsatisfactory results due to systematic deviations from the linear relationship. The deviations can be explained by the formation/decomposition of cementite. Further analysis of the data set shows that the catalyst becomes saturated with carbon after approximately 20 h. In the time interval 25.42–35.09 h (number of measurements N = 20) the results can be approximated by a linear relationship. In this way, it is possible to determine the average reaction rate constant $k_1=4.20·{10}^{-5}\hspace{0.17em}{\mathrm{s}}^{-1},$ $g\left(\mathsf{\alpha }\right)=0.1512\cdot t-2.9254$ ($t$ in hours) ($r^2=0.9891)$. Comparing the values of the kinetic constants determined using Equations (43) and (56), $k_1$ calculated using model F1 stands out as an underapproximation $k_1=2.69·{10}^{-5}\hspace{0.17em}{\mathrm{s}}^{-1},$ $g\left(\mathsf{\alpha }\right)=0.1065\cdot t-2.0956$ ($t$ in hours) $(r^2=0.9882$).

The figure above shows a very small degree of conversion. This corresponds to the maximum distance from equilibrium (DFE), which under extreme conditions leads to the following set of relations:

$$\begin{gathered} DFE=1-z\alpha =0, z=0 \mathrm{a}\mathrm{n}\mathrm{d} \mathrm{t}\mathrm{h}\mathrm{u}\mathrm{s} DEF=1, \\ \mathrm{w}\mathrm{h}\mathrm{i}\mathrm{l}\mathrm{e} \mathrm{f}\mathrm{o}\mathrm{r} \alpha ={\alpha }_{eq}, DFE=0 \end{gathered}$$

(57)

The concept of Distance From Equilibrium (DFE) exists in the literature, but its formal representations vary [117,118]. The basic, definitional expression of DFE is represented by the term given in parentheses in Equation (45), $\left(1-K_x/K\right)$, while Equation (57) provides a justified measure of DFE in the geometric sense.

A separate issue is the degree of conversion achieved close to catalyst deactivation $\left({\alpha }_{\max }\le {\alpha }_{\mathrm{eq}}\right)$. An example of this problem is shown in Figure 5, where the experimental data from the study [38], after recalculation using Equation (51), show the kinetics of the final stage of catalytic activation and a sharp transition to the kinetics of catalytic deactivation. This transition is associated with the inversion of the slope coefficient through the appearance of a change in sign.

It should be noted that in this research the methane stream fed to the reactor was diluted with nitrogen in a 1:1 ratio ($x_{\mathrm{C}{\mathrm{H}}_4\_0}=0.5$), the process was run at $t=750 °\mathrm{C}$ and ${\alpha }_{eq}=1$.

The maximum degree of conversion occurs here around the 120th minute of the test, indicating a transition from the time scale of seconds (for the elementary reaction) to that of minutes or even hours (for the deactivation part). Equation (51) remains correct whenever it takes values higher than the constant term which is determined for the time $t=0$.

From the comparison of the kinetic constants, it can be seen that under flow conditions and with the use of a catalyst, the elementary reaction is unnoticeable on technically usable time scales.

4.4. KCE and EEC

The following considerations are made to assess the extent to which there is a relationship between the parameters of the Arrhenius law and Eyring’s TST. Proof of such a relationship will be beneficial as it will allow comparison of the Arrhenius kinetic parameters determined in different process conditions using the form of the compensation effect KCE, while further transformation of the Arrhenius kinetic parameters into thermodynamic activation functions through EEC will allow analysis of their connections with the second law of thermodynamics.

In order to perform this task, a large set of literature data was used. The parameters of the Arrhenius law were determined from the data given in [13,55,88,102,112,116,119,120,121,122]. The results obtained in this way are presented in the following Figure 6 with a coordinate system of $\left(\ln A\right)$ vs. $\left(E\right)$. This plot includes only those data for which the dimensions of the kinetic constant ($A$) could be assumed to be $\left({\mathrm{s}}^{-1}\right)$. In the process of collecting the data, the pressure unit “atm” was omitted, since it can be related to $\left(P^{\circ }\right)$, while the geometric dimensions were assumed for the given (or estimated) quantities. The error of this approach diminishes when referring to the same conditions, while its overall influence on the precision of the following analysis is marginal due to the use of the logarithmic scale for the frequency factor.

The coordinates of the graph shown in Figure 6 determine a linear relationship of the KCE for the elementary reaction (1). This relationship has a slight asymmetry in the range $E=150-250\hspace{0.17em}\mathrm{kJ}\cdot {\mathrm{mol}}^{-1}$, which corresponds to catalysed reactions. From the top, the line of KCE is limited by the energy of the $\mathrm{C}-\mathrm{H}$ bond. According to [8] its value is equal to $434 \mathrm{kJ}\cdot {\mathrm{mol}}^{-1}$, which corresponds to coordinate (420.7; 36.23). The data were obtained from experiments carried out in a completely mixed reactor with bypass (CPMR) [116]. As a general remark, it is often understood that activation energies $E>250\hspace{0.17em}\mathrm{kJ}\cdot {\mathrm{mol}}^{-1}$ are characteristic of high-temperature pyrolysis processes initiated by the cleavage of $\mathrm{C}-\mathrm{H}$ bonds.

Analogous to KCE, by expressing the Arrhenius parameters in terms of thermodynamic activation functions, it is possible to represent them in terms of EEC. The EEC line for reaction (1) is shown in Figure 7 below [62].

The above expression for EEC is given in the form proposed by [95]:

$$∆S^{+}={\displaystyle \frac{∆H^{+}}{T_{iso}}}+∆S_B$$

(58)

where entropy of activation in isokinetic temperature $∆S_B$ is equal to:

$$∆S_B=\mathrm{R}\ln \left({\displaystyle \frac{k_{iso}}{\mathrm{B}{T}_{iso}}}\right)=\mathrm{R}\ln \left({\displaystyle \frac{A_{E=0}}{\mathrm{B}{T}_{iso}}}\right)$$

(59)

Importantly the components of Equation (58) were determined based on the following relations:

$$∆H^{+}=E-\mathrm{R}{T}_{eq}$$

(60)

$$∆S^{+}=\mathrm{R}\left[\ln \left({\displaystyle \frac{A}{\mathrm{B}{T}_{eq.}}}\right)-1\right]$$

(61)

Equations (60) and (61), when introduced into Equation (58), give KCE with such a difference that for the activation function, the equilibrium temperature $T_{eq}=815.2 \mathrm{K}$ is introduced (Table 1, No. 1). The reason for this is that if $T_{iso}$ is used here, we end up with a tautological problem; therefore, the relationship with the equilibrium temperature was deliberately chosen and is proposed here, since these two values are very different $T_{iso}>T_{eq}$. With this in mind, we obtain:

$$\ln A={\displaystyle \frac{E}{\mathrm{R}{T}_{iso}}}+\ln A_{E=0}+\left[1-{\displaystyle \frac{T_{eq}}{T_{iso}}}+\ln \left({\displaystyle \frac{T_{eq}}{T_{iso}}}\right)\right]$$

(62)

The expression in square brackets in Equation (62) is 0 for $T_{iso}=T_{eq}$ For example, if $T_{iso}=2T_{eq}$, the dimensionless sum is approximately −0.2, and for $T_{iso}>T_{eq}$, the calculated value adds to the intercept. Therefore, the coefficients of KCE and EEC are different because of the reference temperatures used here (intercept) and the regression analysis from the inversion of which $T_{iso}$ is determined (slope).

By substituting in Equation (59), the data given in Figure 6 ($\ln k_{iso}=-2.024)$ and in Figure 7 $T_{iso}=1428.6 \mathrm{K}$, we obtain $∆S_B=-274.80 \mathrm{J}·{(\mathrm{m}\mathrm{o}\mathrm{l}·\mathrm{K})}^{-1}$. Importantly this value is close to the intercept determined in the correlation shown in Figure 7.

Because:

$${\left({\displaystyle \frac{\partial {∆}_r{G}^{\mathrm{\varnothing }}}{\partial T}}\right)}_P=-{∆}_r{S}^{\mathrm{\varnothing }}$$

(63)

than differentiation of Equations (16) and (18) in regard to temperature gives that:

$${∆}_r{S}^{\mathrm{\varnothing }}=∆S^{+}-{∆}_{-1}{S}^{+}$$

(64)

Using Equation (61) twice on the right side of Equation (64), the term $\ln \left(\mathrm{B}{T}_{eq}\right)$ disappears and we obtain:

$${∆}_r{S}^{\mathrm{\varnothing }}=\mathrm{R}\ln \left({\displaystyle \frac{A}{A_{-1}}}\right)$$

(65)

In further considerations, we will assume that $T_{\mathrm{iso}}$ is the same for both “forward” and “reverse” reactions. This is a practical approximation, since we are considering two reactions that take place outside the equilibrium temperature $T_{\mathrm{eq}}$ (where $k$ would be equal to $k_{-1}$). Based on the data given in [62], where the study was carried out under isothermal conditions $\left(T=1073-1313\hspace{0.17em}\mathrm{K}\right)$, for the frequency factor $\left(A\right)$ expressed in $\left({\mathrm{s}}^{-1}\right)$, using Equation (65) we can determine that ${∆}_r{S}^{\varnothing }=100.4 \mathrm{J}\cdot {\left(\mathrm{mol}\cdot \mathrm{K}\right)}^{-1}$. For comparison, according to [1] for $T_{eq}=815.2\hspace{0.17em}\mathrm{K}$, ${\Delta }_r{S}^{\circ }=100.5\hspace{0.17em}\mathrm{J}\cdot {\left(\mathrm{mol}\cdot \mathrm{K}\right)}^{-1},$ while for $T=298.15 \mathrm{K}, { \Delta }_r{S}^{\circ }=80.89\hspace{0.17em}\mathrm{J}\cdot {\left(\mathrm{mol}\cdot \mathrm{K}\right)}^{-1}$. The above results are consistent with expectations when compared with the slope of Equation (20). Within the range of temperatures studied in [62], i.e., from $T_{eq.}=1073 \mathrm{K}$ to $1313 \mathrm{K}$, we observe only slight variability in the increase in entropy with temperature.

Using the reference temperature $T=T_{eq.}$, one can utilize the equality $k=k_{-1}$ which implies:

$$\ln \left({\displaystyle \frac{A}{A_{-1}}}\right)={\displaystyle \frac{E-E_{-1}}{\mathrm{R}{T}_{eq}}}$$

(66)

Combining Equations (65) and (66) with the condition that for ${∆}_r{G}^{\varnothing }=0$, ${∆}_r{H}^{\varnothing }=T_{eq.}·{∆}_r{S}^{\varnothing }$, leads to the following simple relation:

$${∆}_r{H}^{\mathrm{\varnothing }}=E-E_{-1}$$

(67)

For complementarity with any given temperature, Equation (67) can be derived also by combining Equation (16) with Equation (64) to form the relationship:

$${∆}_r{H}^{\mathrm{\varnothing }}=∆H^{+}-{∆}_{-1}{H}^{+}$$

(68)

Using the data given in [62] and Equation (67), we obtain ${\mathsf{\Delta }}_r{H}^{\circ }=124.3\hspace{0.17em}\mathrm{kJ}\cdot {\mathrm{mol}}^{-1}$. The same value (only with a negative sign) is given in work [62] and it is also the intercept in Equation (20). According to [1], for $T_{eq}=815.2\hspace{0.17em}\mathrm{K},$ ${\mathsf{\Delta }}_r{H}^{\circ }=87.48\hspace{0.17em}\mathrm{kJ}\cdot {\mathrm{mol}}^{-1},$ and this value is close to the mean enthalpy given in Table 1, No. 1. The difference between these values is directly due to the inclusion of the solid phase, symbolically assigned to the $\mathrm{F}\mathrm{e}+\mathrm{C}/\mathrm{F}{\mathrm{e}}_3\mathrm{C}$ components, and the inability to distinguish the composition of the solid solution of $\mathrm{C}$ in Fe from the cementite, which is forming and decomposing in undefined stoichiometric ratios.

Comparing Equation (58) with (64) and (67) leads to the interpretation of the thermal and catalytic decomposition of methane from the perspectives of the thermodynamics of activation and phenomenological thermodynamics. Furthermore, Equation (67) is particularly important for the reversible reaction because it satisfies the empirical observation $(E>{\mathsf{\Delta }}_r{H}^{\circ })$.

In EEC, the intercept, which according to Equation (59) is the constant of the activation entropy and takes values $<0$, occurs when $k_{\mathrm{iso}}=A_{\left(E=0\right)}<\mathrm{B}{T}_{\mathrm{iso}}$. According to Equation (64), the difference between the activation entropy and the initial state entropy follows directly from the inequality $\mathsf{\Delta }{S}^{+}>{\mathsf{\Delta }}_{-1}{S}^{+}$ and because ${\mathsf{\Delta }}_r{S}^{\circ }>0$. Furthermore, according to Equation (67) $E>E_{-1}$.

Comparing Figure 6 and Figure 7 we can conclude that we are dealing with the same information. For both relationships $r^2=0.898$ for N = 35 measurements. Furthermore, both are significant with a probability of 99.9%. The fundamental difference between the two lies in their interpretation. From the determined EEC line and Equation (58), the increase in the activation entropy in the range $\mathsf{\Delta }{S}^{+}=-200-100\hspace{0.17em}\mathrm{J}\cdot {\left(\mathrm{mol}\cdot \mathrm{K}\right)}^{-1}$ can be explained according to Equation (61) by the increase in the reaction rate, which according to Arrhenius leads to a higher pre-exponential constant, i.e., $(\mathsf{\Delta }{S}^{+}\propto \ln A)$.

This may also explain the red point below the EEC line in Figure 7 (not included in the correlation). It can therefore be argued that it indicates a significant decrease in the rate of methane decomposition under conditions of carbon absorption/dissolution by the excess solid phase $(\mathsf{\gamma }-\mathrm{Fe})$. The conclusion is that while we use correlations expressed in terms of ratios, here Equations (18) and (65), the problem is different for absolute expressions, e.g., here Equations (16) and (64). Equation (19) is an example of this, because the decrease in methane decomposition kinetics is related to the course of the reaction or physical process: $(3\mathrm{Fe}+\mathrm{C}\leftrightarrow {\mathrm{Fe}}_3\mathrm{C})$ (Table 1, No. 5), which only translates the formation of weaker bonds in the solid solution.

Looking at the reaction path (19) in the opposite direction, we have no guarantee that we will return to the structure symbolically written as $\left[\mathrm{C}\right(\mathrm{i}\mathrm{n} \mathsf{\gamma }-\mathrm{F}\mathrm{e}\left)\right]$, since $(\mathsf{\gamma }-\mathrm{Fe})$ cannot be a product of reaction (1).

The occurrence of a physical transformation in the system of the process considered here affects its kinetics in the context of the possibility of analysing its reversibility. An illustration of the Bell–Boudouard reaction (Table 1, No. 12) can be found in Figure 1 in [123], since the closing of the chemical cycle does not always lead to a return to identical structural forms.

In line with the authors’ previous work [95], such an analysis leads to the following logical sequence of relationships: kinetics → activation thermodynamics → phenomenological thermodynamics. Importantly, this provides an analytical link between kinetics and thermodynamics.

5. Discussion

Now that the difficulties of modelling the process of methane pyrolysis under different initial process conditions are clear, the discussion of the current understanding of thermodynamic and kinetic modelling can finally be presented in a more holistic manner. As mentioned above, it is useful to distinguish between two different stages of the process, i.e., activation and deactivation. For the reasons given in the introduction to this review, the following discussion is based on the assumption that reaction (1) is catalysed by iron in the form of iron ore, with haematite as the reference material.

5.1. Range $\alpha \to {\alpha }_{eq}$

In the initial moment of the reaction $\alpha =0$ while $DFE=1$. On the other hand, at its end $\alpha ={\alpha }_{\mathrm{eq}}=1$ and $DFE=0$. It has already been shown that the time dependence of α is a linear relationship, but only when the reaction is carried out in the range of small degrees of conversion. In such a case, Equations (43) and (56) are used, but this model cannot be said to describe reaction (1) completely and must therefore be regarded as an approximation. These equations are further simplified when ${\alpha }_{\mathrm{eq}}\to 1$ and such process conditions can be achieved by diluting the methane with an inert gas or by increasing the temperature. Both allow the term responsible for the reversible reaction $k_{-1}=0$ to be ignored.

To further increase the accuracy with which the model based on the above Equations (43) and (56) predicts the course of reaction (1), it must be fitted with additional terms. Such attempts are presented above in Table 4 and Figure 4. Model F1 (No. 5) predicts a significantly different value of the kinetic constant $k_1$ compared to the results obtained with other models presented in the table. This is probably due to the fact that the remaining models are richer in terms that more accurately describe the nature of the process involved in the pyrolysis of methane. By interpreting the differential equations, it can be shown that they approach quasi-second order. For small $\mathsf{\alpha }$, this becomes particularly clear for model No. 2 $\left({\displaystyle \frac{1-\mathsf{\alpha }}{1+\mathsf{\alpha }}}=(1-\mathsf{\alpha })(1+\mathsf{\alpha }{)}^{-1}\cong (1-\mathsf{\alpha }{)}^2\right)$. Further analysis of these models leads to the conclusion that their basis is still the F1 kinetics but modified by terms characteristic of the uncatalysed reaction (1) and phenomena occurring during the transition from homogeneous to heterogeneous reaction. A detailed analysis of the kinetics of this type of reaction is given in [111]. For flow reactors (e.g., plug flow reactor), the reaction time, measured as the residence time, is the ratio of the reactor volume to the flow rate, which is proportional to the integral $t\propto \int {\displaystyle \frac{d\alpha }{r}}$.

Despite the catalysed nature of reaction (1), it is clear that the elementary reaction for the thermal process is still an important element to be considered. Referring to the results shown in Figure 3, the following Arrhenius parameters were presented there: $E=130.5\hspace{0.17em}\mathrm{kJ}\cdot {\mathrm{mol}}^{-1}, \ln A=8.59$, where $A$ is given in $\left({\mathrm{s}}^{-1}\right)$ and $P=5.61\hspace{0.17em}\mathrm{M}\mathrm{P}\mathrm{a}$. These values can be compared with results from the opposite end of the spectrum, e.g., data obtained for higher temperatures $>1200°\mathrm{C}$ and at atmospheric pressure [112,113,114,121], where the following kinetic parameters were determined: $E=420.7\hspace{0.17em}\mathrm{kJ}\cdot {\mathrm{mol}}^{-1}, \ln A=36.23$, where $A$ is given in $\left({\mathrm{s}}^{-1}\right)$. It is important to note that both of the above examples indicate that the reaction takes place on a second or millisecond scale.

According to Equations (8) and (12), in the case of high pressure $P=5.61\hspace{0.17em}\mathrm{MPa}$, at $t=900 °\mathrm{C}$ we obtain ${\alpha }_{\mathrm{eq}}=0.4336$. Thus, the real experimental conversion must be lower than $\alpha <0.43$ and therefore the characteristic Arrhenius parameters determined above must be considered as underestimated.

Another case for comparison can be found in works [116,124], where the former study presents results obtained for three types of reactors and a very wide range of activation energies, while the latter presents the results of computer simulations carried out for $T_{\mathrm{iso}}=1190.9\hspace{0.17em}\mathrm{K}$. It is noteworthy that the rate constants in both articles are the same and therefore we observe the unambiguity of the kinetic constants for all the dependencies related to this temperature.

5.2. Range $\alpha \to 0$

The part of the course of methane pyrolysis that takes place under conditions of decreasing conversion is sometimes referred to as the deactivation phase of the catalyst test. Figure 5 above shows this period, but its formal analysis boils down to calculating the amount of carbon deposit formed. As shown above, this information can be determined using Equations (25) and (28). Analysing this figure, we can see that at the maximum (i.e., for $t^{\mathrm{*}}=119.80$ min) according to Equations (8) and (53) ${\alpha }_{max}=0.7864 \mathrm{w}\mathrm{h}\mathrm{i}\mathrm{l}\mathrm{e} {\alpha }_{eq}=0.8211$.

5.3. Range $\alpha \to {\alpha }_{eq}\to 0$

A separate issue that needs attention is whether it is possible to model the kinetics of reaction (1) in the full range of variability of its degree of conversion. To date, no such analyses have been found in the literature, although the proposals presented in [125] may provide some guidance. The simplest form capable of describing the full range of variability of the degree of conversion can be taken as an empirical formula:

$$\alpha ={\displaystyle \frac{k_1}{k_1-k_3}}\left(e^{-k_{3}t}-e^{-k_{1}t}\right), k_1>k_3$$

(69)

The kinetic constants for Equation (69) are presented in [125].

6. Conclusions

A literature review on the kinetics of the thermal decomposition of methane was carried out, focusing on finding the elementary reaction as defined in the IUPAC Gold Book [107], and considering the importance of the reversible term (see the initial approach of Equation (32)). It was found that the limitations on the efficiency of the process of thermal decomposition of methane are the equilibrium conversion degree ($\alpha$) and the distance from equilibrium (DFE). These limitations can be mitigated by increasing the process temperature to at least 900 °C and diluting the methane with an inert gas. These measures in turn make reaction (1) irreversible.

Detailed analyses of the non-catalytic elementary reaction of the purely thermal process indicate that at high pressures or temperatures and at short process times (on the order of ms), the reaction kinetics can be described by the parameters of the Arrhenius law. Despite the wide range of activation energies $E=130.5–420.7 \mathrm{k}\mathrm{J}·\mathrm{m}\mathrm{o}{\mathrm{l}}^{-1}$ obtained by the calculations presented here, these parameters satisfy the Kinetic Compensation Effect (KCE) and the Enthalpy–Entropy Compensation (EEC) for thermodynamic activation functions. It is noteworthy that many descriptions of catalytic processes found in the literature indicate a much longer time scale (with durations of the order of h). However, these experiments cannot be considered to have been carried out under conditions that are well adapted to the characteristics of the reaction presented here (1).

It has been shown that in the presence of catalysts, the kinetic equations covering the range $\alpha \to {\alpha }_{max}\to {\alpha }_{eq}$ can take simple forms, such as the F1 model. After supplementing this approach with terms related to the characteristics of the process (effect of purely thermal conversion, reversibility, heterogeneous reaction, etc.), this model exhibits a quasi-second order (Nos. 1–4 in Table 4). Furthermore, $\alpha$ is also a very convenient kinetic parameter, particularly useful for balancing the process and determining the amount of carbon deposit produced.

The potential for utilising the Kinetic Compensation Effect (KCE) to elucidate the underlying mechanisms of reaction (1) and to discern both synergistic and antagonistic catalytic interactions occurring within this process was demonstrated. The phenomenon appears to manifest itself following the transformation of the Arrhenius kinetic data into thermodynamic activation functions (EEC). It is important to note that there are a number of potential reasons for the existence of compensation relationships. In some cases, the underlying cause may be a fundamental relationship, while in other instances, the phenomenon may be attributed to experimental errors or technical limitations. Moreover, Figure 7 offers preliminary indications of a potential subsequent relationship between the analysed literature data. However, further proof can only be provided once a larger sample size is available.

The thermal decomposition reaction occurring in the system containing: methane, carbon, and iron at temperatures close to the isoconversion temperature depends on the proportions of the components and the effects of the interactions between the Fe/C phases. Elements of Transition State Theory (TST) were applied here to the analysis of the system by comparing the ratio of the “forward-backward” kinetic constants with the phenomenological thermodynamic functions. Here, quantities of phenomenological thermodynamics, i.e., Gibbs free energy, enthalpy, entropy (${∆}_{\mathrm{r}}{\mathrm{G}}^{\mathrm{\varnothing }},{{∆}_{\mathrm{r}}\mathrm{H}}^{\mathrm{\varnothing }},{{∆}_{\mathrm{r}}\mathrm{S}}^{\mathrm{\varnothing }}),$ were directly related to activation functions ($∆{\mathrm{G}}^{+},∆{\mathrm{H}}^{+},∆{\mathrm{S}}^{+})$.

Finally, for iron-based supported catalysts, it is necessary to consider their activation. The stages of catalyst preparation, apart from the necessary drying, grinding, or sieving, are more important from a process point of view in choosing the conditions under which the subsequent reduction takes place (e.g., with methane or hydrogen). In the context of methane pyrolysis, the preferred form of iron catalyst is $\mathsf{\alpha }-\mathrm{F}\mathrm{e}$, which is also the most favourable structure of Fe for cementite formation. It is important to note that the polymorphic form of $\mathsf{\gamma }-\mathrm{F}\mathrm{e}$, which occurs at temperatures > 912 °C has a lower catalytic activity and therefore its formation should be avoided. It should also be noted that the conversion of the supported iron oxide catalysts, and especially the naturally occurring ores, to their active forms of iron is a time-consuming and reductant-intensive process.

The analysis presented here suggests two main areas for future research. Firstly, it would be beneficial to improve the description of reaction (1) by incorporating Equation (69). This has not been performed previously, but it allows the thermal decomposition of methane to be described in terms of the kinetic constants associated with the thermodynamic activation functions (TST). This may facilitate the integration of current knowledge of the reaction pathway and catalyst deactivation. Secondly, it can be seen that further research is needed on the irreversible [13], dissociative adsorption [110], and irreversible dissociative adsorption [126] models. These should be analysed in the context of reactions (2) and (3), or even described in the light of different conversion pathways, such as acetylene–benzene [46]. Furthermore, fractional order models have not been analysed here, although they may be appropriate. In regard to catalysts other than iron, the literature so far presents a number of cases where activated carbons show significant catalytic activity, but there is no research published yet where classical allotropic forms of carbon are used. The boundary between the carbon-catalysed and the purely thermal process of methane decomposition also needs to be better described. For nickel, it is known that carbides similar to Fe can be found: Ni₃C, Ni₈C [13,45], or Ni₈C on TiC substrate [47], but their analysis in the context of reaction (1) is not yet as well understood. The importance of the reversible term of reaction (1) also needs to be further developed in the context of its impact on the activity of the catalysts towards the carbon deposit formed, as it can be difficult to determine such a model that would fully grasp the phase transformations of the catalyst and its reactions with both the gas atmosphere and the $C$ deposit.