Modelling of the catalytic initiation of methane coupling under non-oxidative conditions

The experimentally observed interplay between catalytic activation of methane on Fe © SiO 2 and gas-phase free radical methane coupling under non-oxidative conditions is analyzed by mechanistic modeling as well as by experiments. For the modeling, an off-the shelf gas-phase model, AramcoMech 3.0, was used unaltered to keep the number of adjustable parameters as low as possible. It was complemented by surface reactions specifically accounting for methane activation to methyl radicals. The model was validated against an independent set of experimental data and exhibited good accordance. The model accurately captured the significant contribution of gas-phase reactions responsible for methane conversion in the post-catalytic zone, indicative of gas-phase autocatalytic methane coupling. The low-activity induction period in gas-phase methane pyrolysis can effectively be overcome by adequate catalytic activation. Results show that the catalytic reaction only influences the activity of the system, with gas-phase reactions dictating the selectivity distribution. Simulations demonstrated that the optimum catalytic conversion roughly amounts to 4 % at 1000 ◦ C and 1 atm. An equivalent effect can be reached by adding ca. 2 % of ethane or ethylene to the feed. Detailed reaction-path analyses were employed to corroborate these phenomena. Gas-phase reactions were found to be very rapid at 1000 ◦ C, hence determining the product selectivity, without impact from either catalyst or C 2 hydrocarbon addition. Current, freely available gas-phase models lack the required accuracy for detailed kinetic modeling of the product distribution, showing the requirement for the development of a dedicate non-oxidative methane coupling model.


Introduction
Natural gas, consisting of 75-99 % of methane depending on its origin [1], is considered a strategically useful source for C 2+ hydrocarbon production [2].The large quantity of proven natural gas reserves [3] combined with the expected decrease in crude oil production [4] leads to the prediction that natural gas will supersede crude oil as the main source of olefins and aromatics in the 21st century [5].Industrial processes are already operated at large scale for the synthesis of fuels and base chemicals, e.g., by the Fischer-Tropsch [6] and methanol-toolefins [7] reactions.These processes all involve multiple conversion steps, starting with methane steam reforming for synthesis gas (CO and H 2 ) production [8] and are only viable at large installed capacities [9].Direct methane coupling to higher hydrocarbons has therefore been a highly researched topic in the chemical community for over a century [10,11].
Early research focused on non-catalytic conversion of methane at high temperatures [10,12], i.e. 'methane pyrolysis'.The main drawbacks of methane pyrolysis are these high operation temperatures (950-1700 • C) and/or long residence times (6-15 s @ 1000 • C) required to achieve significant methane conversion [10,13].These challenges stem from the high thermal stability of methane [14] as well as the pronounced endothermicity of the dehydrogenative coupling reactions [15].Other important issues in methane pyrolysis are significant carbon deposition at high methane conversion, high energy consumption because of the high operation temperature as well as relatively low product yields, implying difficult separation [16][17][18].Starting in the 1970 s, a tremendous effort was made to understand the highly complex free radical reactions involved in methane pyrolysis.Chen et al. [19,20] reported on the reactions responsible for the formation of the primary product, ethane, and secondary product, ethylene.Their model showed great compliance with experimental data at low methane conversion.Roscoe and Thompson continued this investigation [21] and were the first to propose a mechanism explaining the autocatalysis phenomenon, in which the methane conversion rate increases sharply after a so-called induction period.These authors attributed this phenomenon to the coupling reaction involving ethylene, acetylene and hydrocarbons with more than two C-atoms.This model was further developed by Dean [22], who showed the importance of the reversibility of the reactions proposed by Roscoe and Thompson [21].This work was continued in the 90 s by the research group of Holmen [15,[23][24][25] as well by Billaud and Guéret [14,[26][27][28][29].The culmination of this work has been summarized in a model by Matheu et al. [30], which gave extra attention to the contributions of propyne, allene and fulvene in the autocatalytic cycle.Oxidative coupling of methane (OCM) [31,32] as well as methane dehydroaromatization (MDA) [33,34] are alternative strategies to overcome the challenges encountered in methane pyrolysis, but both still result in too low of a product yield for industrial application.
In 2014, Guo et al. [35] reported a Fe/SiO 2 catalyst capable of coupling methane non-oxidatively to ethylene and aromatics at temperatures above 950 • C.This catalyst is referred to as Fe©SiO 2 to denote the lattice confined single atom nature of the iron.Most crucially, these authors reported the absence of coke formation, even at high methane conversion (48 %), and catalyst stability up to 60 h.These performance figures proved difficult to reproduce [36][37][38][39][40][41], nevertheless these latter publications contributed to better understanding of the catalytic system.The sole function of the Fe©SiO 2 catalyst in the overall reaction mechanism is the formation of CH 3 • radicals, while all other reactions, including the coupling reaction, occur in the gas-phase [35].Another DFT study by Kim et al. [39] confirms this assessment, although these authors propose that direct acetylene formation on the catalyst is theoretically also possible.Liu et al. [42] propose that C-C coupling to directly form ethylene can occur on the Fe©SiO 2 catalyst.The work by Toraman et al. [43], for the first time, combined a detailed surface mechanism with an extensive propene pyrolysis gas-phase model from Wang et al. [44] to elucidate the interplay between heterogeneously catalyzed and gas-phase reactions.These authors propose that the Fe©SiO 2 catalyst predominantly forms ethylene, with catalytic CH 3 • radical generation playing only a minor role.Thus far, no experimental evidence has been provided to prove either the catalytic CH 3 • or C 2 H 4 formation hypothesis.
The significance of gas-phase reactions was experimentally validated in our previous work [41,45] as well as in a patent filed by Sabic [46].It was shown that the catalyst is required for initiation of the coupling reaction, after which post-catalytic residence time is the predominant factor in propagating the coupling reaction, leading to significant methane conversion in gas-phase.A long residence time inside the catalyst bed and residence time at reaction temperature upstream of the catalyst bed were both detrimental to the total hydrocarbon selectivity, enhancing formation of deposits, without significantly increasing methane conversion.The addition of ethane and ethylene can significantly increase the methane conversion rate [45,47,48].
The present work provides a semi-quantitative assessment of the interplay between the initiation of the reaction on the Fe©SiO 2 catalyst and the subsequent gas-phase coupling reactions, focusing on catalytic CH 3 • radical generation via an in-house developed microkinetic model in combination with an existing gas-phase microkinetic model.Emphasis is put on the impact of the catalytic reactions on the rate of gas-phase reactions, involving a detailed quantification of catalytic initiation and gas-phase propagation of the NOCM reaction.Reactionpath analysis is used to explain the experimental product distribution.The model is used to understand how to maximize reactor productivity by optimizing the catalyst-free volume and catalyst contact time.Lastly, the effect of addition of ethane and ethylene on the methane coupling reaction is elucidated.

Catalyst synthesis
The Fe©SiO 2 catalyst is synthesized according to the method described by Guo et al. [35].The Fe 2 SiO 4 precursor was synthesized according to the method described by DeAngelis et al. [49], which is then mixed with high purity quartz and milled overnight in a planetary ball mill.The resulting powder is fused for 6 h in air at 1700 • C. The resulting slab is crushed and sieved to the desired particle size and leached for 2 h in 0.5 M HNO 3 .Any further details as well as catalyst characterization has been reported in our previous work [41].Two particle size fractions are used, namely 100-250 µm and 250-500 µm.The data obtained with the 100-250 µm fraction is used for fitting of the model because it meets all criteria for intrinsic kinetics measurements following the method described by Berger [50], i.e. mass and heat transport are not limiting and ideal plug-flow behavior is obtained.The 250-500 µm fraction is used for validation and further analysis of the gas-phase-catalyst interaction, because this is the fraction used in previous publications [41,45].Note that the 250-500 µm fraction is predicted to have a slight radial temperature profile, albeit below 3 K.

Kinetic measurements
Catalytic tests are carried out in a quartz tube, ID 4 mm, which is placed in a 3-zone oven with thermal insulation between each zone, described in detail in our previous publication [41].The catalyst is inserted into the desired position in the quartz reactor, as displayed in Fig. 1a.The catalyst is held in place by small quartz wool plugs.The catalyst is first exposed to 90 vol% CH 4 in N 2 at 1000 ml⋅g cat -1 ⋅h − 1 for two hours at 900 • C for activation purposes, following the procedure of Bao pers.commun.[51], to carburize the active site forming a FeC 2 Si site [35].After activation, the temperature is increased to the desired reaction temperature, in the range between 950 • C and 1100 • C. Space velocity is defined in ml⋅g cat -1 ⋅h − 1 , with the volumetric flow taken at standard conditions (25 • C and 1 atm), referring to the combined feedgas flowrate (CH 4 + N 2 + C 2 ).A blank measurement following the same procedure as catalytic measurements was performed to ensure that the reactor is exposed to the same pretreatment prior to measuring.Blank measurements are performed in an empty reactor tube, on similar note, the free volume presented in Fig. 1a refers to gaseous volume.No beddiluent is used due to the activity of any surface area at reaction conditions [16,38].Product mixtures are analyzed using a three-channel Varian CP3800 online gas chromatograph (GC).Any further details concerning the setup, experimental procedures and calculation methods can be found in our previous work [41].The temperature profile considered in the model is based on the measured temperature profile in the lab-setup (SI Fig. S1).

Model development
A top-down approach has been chosen for the model development due to the system's complexity.Firstly, the physical description of the reactor is generated, due to the significant impact of the temperature profile on the reaction [41].Afterwards the gas-phase model is selected based on its compliance with the blank experiment and lastly the catalytic surface reactions are included.

Reactor model
ChemKin-Pro 2020 R2 was used to perform the calculations.The constructed reactor model comprises four consecutive zones, simulated as plug-flow-reactors (PFR), see Fig. 1b.The first and last zones model the pre-heater and post-heater respectively and exclusively consider gas-phase reactions.Zones 2 and 3 simulate the reactor-zone, by splitting it into a part with catalyst where surface reactions contribute, and a part without catalyst considering exclusively gas-phase reactions.The different catalyst positions shown in Fig. 1a can be accurately described in this way.The temperature profile presented in SI (Fig. S1) was used as the temperature profile in the four zones, in order to account for reactions occurring during heating and quenching of the gas stream.The assumption of plug-flow behavior has been validated by comparing the diffusion rate of H⋅ radicals, calculated via the Fuller, Schettler and Giddings's method [52] with the reactor characteristics, according to the criteria presented in [50].The analysis presented in supporting information section SI. 2, show that plug flow can be assumed for all cases, except the most extreme, i.e. the longest residence time and the largest particle size fraction.

Gas-phase model
Six gas-phase models, available in literature, i.e.Matheu et al. [30], Aramcomech 3.0 [53], Nuigmech 1.1 [54], Gudiyella et al. [55], Chu et al. [56] and Wang et al. [44] were evaluated based on their compliance with the blank experiments (see supporting information S.3).It must however be noted that none of the tested gas-phase models gave a quantitatively satisfactory prediction over the measured range of methane conversion.It is thus suggested that the current models should either be updated to accurately predict methane pyrolysis or a specific methane pyrolysis mechanism should be developed.Still, AramcoMech 3.0 [53] provided an accurate description of methane autocatalysis, the accurate prediction of the product range as well as reasonable trends in methane conversion and product selectivity overall and was chosen.Though it must be stated that the models of Gudiyella et al. [55], Chu et al. [56] and Wang et al. [44] show the same autocatalytic mechanism and could also have been chosen as alternative to AramcoMech 3.0.A detailed comparison of the mechanisms is provided in supporting information S.4.We found that the choice of gas-phase model between these four has no impact on the main trends and thus on the conclusions thereof.
The two major shortcomings of the AramcoMech 3.0 model in reproducing the experimental observations are situated in the ethylene selectivity, which is significantly overpredicted, and in the selectivity towards the 2-ring aromatics, naphthalene and indene.The overprediction of the ethylene is characteristic of the AramcoMech 3.0 model.Concerning the 2-ring aromatics, the model seems to replace the naphthalene (C 10 H 8 ) selectivity with indene (C 9 H 8 ), which is actually not detected at all.To address this discrepancy, it is assumed that indene predicted by the model represents naphthalene observed in the measurements, including a correction for the difference in number of carbon atoms in the molecules.The parity diagrams presented in SI Fig. S10 show that the gas-phase kinetic model predicts the selectivity to the major products well.

Surface reaction model
The catalytic surface reactions are directly adopted from the DFT calculation shown in the work of Guo et al. [35] and the catalytic cycle is schematically presented in SI Fig. S9.In short, the cycle involves dissociation of two methane molecules to methyl radicals, which desorb, and two adsorbed H radicals, which couple to form molecular hydrogen.The activation energies are based on the energy barriers presented in the original publication [35] and the pre-exponential factors are manually adjusted and evaluated to match experimental results.More details concerning the mechanism as well as the adjustment method can be found in Supporting information section S.5.

Model validation
The pre-exponential factors of the rate coefficients of the catalytic surface reactions resulting from the parameter adjustment to match the 100-250 µm fraction measurements are shown in supporting information Table S1.Fig. 2 shows the parity plots comparing the developed model against the experimental data, whereas Fig. S11 in the supporting information presents the same data in linear plots.It is clear from Fig. 2a and c that the model accurately predicts the methane conversion when varying either the residence time (R 2 = 0.992) at constant temperature, or the temperature at constant residence time (R 2 = 0.999).The deviation at higher conversion can be attributed to the absence of coke formation in the model, which in reality will partly block active sites, thus limiting conversion.It must however be noted that gas-phase freeradical chemistry dominates and the agreement between the experimental and modeled methane conversion is only determined by the level of methane activation by the catalyst.Therefore, it is impossible to experimentally distinguish between the two DFT based catalytic mechanisms proposed by Guo et al. [35] and Toraman et al. [43] respectively.The product selectivity presented in Fig. 2b and d also shows reasonable agreement, although the selectivity to renowned coke precursor molecules such as acetylene, benzene and naphthalene are overpredicted due to the absence of coking in the model.Ethylene selectivity is also Validation of the model was performed against experimental data on the Fe©SiO 2 catalyst using a larger particle size fraction, i.e. 250-500 μm [41,45].The active surface area per length of the catalyst bed was thus reduced in the model from 129 mm 2 ⋅mm − 1 to 60 mm 2 ⋅mm − 1 , accounting for the difference in average particle size in the two data sets.The results of the validation are shown in Fig. 3 (a copy of these graphs without logarithmic axis can be found in Fig. S12).Clearly, the experiments with varying space velocity and temperature (Fig. 3a) are also accurately described for the larger particle size fraction (R 2 = 0.990).The product selectivity is predicted even better for the larger particle fraction, due to the lower coke selectivity observed with the higher particle size fraction [41].

Catalytic initiation of the reaction
Gas-phase reactions contribute significantly to both product formation as well as methane conversion [41,43].To understand their role, simulations are performed evaluating the effect of free volume at 1000 • C upstream and downstream of the catalyst bed, accounting for both gas-phase reactions and catalytic reactions and compared with the corresponding experiments.The effect of free volume upstream or downstream of the catalyst (Fig. 1) is accurately predicted by the model, as shown in Fig. 4a.Product selectivity is also predicted reasonably well as shown in Fig. 4c and d.Fig. 4b shows the modelled effect of total gasflowrate on methane conversion with free volume upstream or downstream of the catalyst bed, using the full reactor model.The inverse of the flowrate on the x-axis is proportional to both residence time inside the catalyst bed as well as the residence time in the gas-phase upstream or downstream of the catalyst.At high flowrate, i.e. short residence times in both free volume and catalyst bed, there is no difference Fig. 2. Parity plots comparing the results of the model with the experimental results measured with the reactor completely filled with 840 mg Fe©S-iO 2 100-250 µm (Fig. 1a); pre-heater and post-heater at 400 • C; (a) Methane conversion and (b) product selectivity to major hydrocarbons when changing the SV between 240 and 4000 ml/g cat /h with the reactorzone at 1000 • C; (c) Methane conversion and (d) product selectivity to major hydrocarbons when changing the temperature of the reactor zone between 950 and 1100 • C, at constant SV of 2400 ml/g cat /h.Fig. 3. Parity plots comparing the results of the model with the experimental results measured with the complete reactor zone (Fig. 1b) filled with 840 mg Fe©SiO 2 250-500 µm; pre-heater and post-heater at 400 • C (Fig. 1a), previously published in [41] (a) Methane conversion and (b) product selectivity to major hydrocarbons when changing the SV between 240 and 4000 ml/g cat /h with the reactor-zone at 1000 • C and when changing the temperature of the reactor zone between 950 and 1100 • C, at constant SV of 2400 ml/g cat /h.  between extended gas-phase residence time upstream (B2) or downstream (T2) of the catalyst bed.At lower flowrates, however, the conversion with extended gas-phase residence time downstream of the catalyst bed (T2) becomes significantly higher compared to extended gas-phase residence time upstream (B2).These results are in agreement with experimental results in our earlier publication as well [41].The significant increase in conversion in gas-phase in the post-catalytic volume is the result of catalytic initiation of the free-radical coupling reaction.This is in agreement with the findings of Toraman et al. [43] who showed that methane conversion increased from 6 % directly after the catalyst bed to 13 % in the post-catalytic volume in gas-phase model they employed.

Contribution of methane conversion in gas-phase
To simplify the analysis and better understand the catalyst's role, an idealized system, i.e. assuming an uniform reactor temperature of 1000 • C and no conversion in the volume upstream of the catalyst, was assessed.Fig. 5a shows the effect of residence time (analogous to 1/ flowrate in Fig. 4b) on methane conversion, equivalent to an axial concentration profile, for the T2 system in which the catalyst contributes via catalytic initiation.Fig. 5b shows the conversion curve in a purely gas-phase system, i.e. without catalytic initiation.The T2 system (Fig. 5a) exhibits significantly more conversion in the free-volume section downstream than the empty reactor (Fig. 5b) at the same residence time (note the different scales of the time axis in Fig. 5a and 5b).
In order to assess the impact of the catalyst in the gas-phase reactions, the conversion at the end of the catalyst bed (about 4 %, indicated with the dotted line) is matched with the same conversion in the blank experiment (in practice, by combining graphs in Fig. 5a and b in Fig. 5c, with the x-axis of Fig. 5a shifted to the right).The same method is used for all graphs in Fig. 5c and d as well as Fig. 6.Fig. 5c and d and Fig. 6 thus show the conversion and concentration profiles as a function of residence time for the T2 case as well as the empty reactor.
The simulation of the empty reactor shown in Fig. 5b shows a gradual shift in regime: methane conversion rate at short residence time is low but increases with increasing residence time and conversion, typical for slow induction followed by autocatalysis.Fig. 5c shows that the catalyst effectively shortens the induction phase, while not influencing the conversion profile in the post-catalytic free volume.This was further verified by varying the catalyst residence time (see Fig. 5d).Note that this would also be the case if the conversion in the catalyst bed would be increased by using a catalyst more active in methyl radical formation.Fig. 6d also indicates that the influence of the catalyst on conversion is the largest when aiming at roughly 4 % conversion.This is analyzed in more detail mathematically in supporting information S.8.Lower conversion in the catalyst bed would lead to a significantly larger reactor because more volume downstream of the catalyst would be required, whereas the marginal gain in reactor volume becomes very small with higher conversion in the catalyst bed.However, the main disadvantage of a larger catalyst bed is increased formation of carbonaceous deposits [41].

Mechanistic evaluation of acceleration in methane conversion rate
To understand the catalyst's role on the kinetics, Fig. 6 shows the concentration profiles of important intermediate species and products as a function of the residence time, similar to Fig. 5 for methane conversion.C 1 , C 2 and H⋅ species are the main contributors to the chemistry involved in gas-phase methane activation [21,22].The concentration Fig. 6.Concentration profiles of (a) ch 3 ⋅ radicals (b) C 2 hydrocarbon (c) C 2 H 5 ⋅ and C 2 H 3 ⋅ radicals (d) H⋅ radical; concerning the catalyst-initiated system and a pure gas-phase system at a constant temperature of 1000 • C as in Fig. 5.

R.S. Postma et al.
profiles as a function of the residence time in the empty reactor as shown in Fig. 6 are based on the microkinetic model for gas-phase methane auto-catalysis according to the AramcoMech 3.0 model [53].The schematic overview of this reaction cycle can be found in Fig. 7, this is an extended version of the cycle presented in our previous work [45].The following section details this mechanism with respect to the concentration profiles for the empty reactor in Fig. 6.
Firstly, the profiles in Fig. 6 will be discussed in light of gas-phase chemistry mentioned in earlier publications [19][20][21][22] (see reaction (1) to reaction (10) below), followed by a discussion on the interplay between catalytic and gas-phase reactions.In the induction phase, i.e. at low methane conversion and in absence of catalysts, slow activation of methane in the gas-phase (reaction (1)) dominates.It cannot be ruled out that the reactor-wall or impurities are involved.The potential heterogeneous nature is included in the model in reaction ( 1) by (m) active site species.Initially each H⋅ radical formed via reaction (1) will activate another methane molecule via reaction (2), causing equal rates of reactions ( 1) and ( 2) during the induction phase.The H⋅ radical concentration stays low even though that of CH 3 ⋅ increases significantly, see Fig. 6a and d, due to the rapid consumption of H⋅ radicals in reaction (2).The CH 3 ⋅ radicals combine to form ethane (reaction (3)); subsequently, ethane reacts with another CH 3 ⋅ radical, forming an C 2 H 5 ⋅ radical and CH 4 (react.( 4)).This reaction is fast, as demonstrated by the low ethane concentration in Fig. 6b, while both C 2 H 5 ⋅ radical and ethylene concentration keep increasing (Fig. 6b and c).C 2 H 5 ⋅ radicals are very reactive and readily dissociate to ethylene and a H⋅ radical (eq.( 5)), explaining the rapid and immediate increase in ethylene seen in simulations in Fig. 6b in absence of catalyst.Similar to ethane, ethylene can dehydrogenate further towards acetylene following reactions ( 6) and ( 7) with C 2 H 3 ⋅ as reactive intermediary radical.Note that, even though C 2 dehydrogenation causes the formation of H⋅ radicals, the consecutive dehydrogenation of C 2 hydrocarbons, shown in eqs.( 4)-( 7), has no effect on methane conversion, since for every produced H⋅ radical (reactions ( 5) and ( 7) another methane molecule is generated from a CH 3 ⋅ radical (reactions (4) and ( 6)).This shows the need for carbon-carbon coupling to accelerate the auto-catalytic cycle.Ethylene is significantly less reactive than ethane and CH 3 ⋅ and C 2 H 5 ⋅ radicals as is apparent from Fig. 6b, showing significant higher concentrations of ethylene and acetylene compared to ethane.The Aramcomech3.0 significantly overpredicts the ethylene concentration, as noted in section methodology section on model development, as well as shown in SI sections S. 3 and S. 10.This will likely lead to a mis-match between the model results and the experimental results in the acceleration of the autocatalytic cycle, giving an explanation for the earlier onset of the autocatalytic acceleration in the measured data as shown in SI Fig. S3.This discrepancy will diminish when the autocatalytic cycle is fully established and the formation and consumption rate of ethylene balance.
Note that, even though H⋅ radicals do not selectively react with methane, the reaction with methane is dominant because of the high methane concentration compared with all other molecules and radicals.For the autocatalysis to take effect, a net formation of H⋅ radicals needs to occur [21,22], to activate methane via reaction (2).All reactions generating these surplus H⋅ radicals require C -C coupling to higher hydrocarbons.Reaction (1) is an exception, but this reaction is slow and constant in gas-phase in absence of catalyst.Acetylene exclusively reacts to higher hydrocarbons according to the model, whereas in reality it is also subject to dehydrogenation resulting in coke formation [57].The formed C 3+ hydrocarbons are very reactive under reaction conditions and spontaneously dehydrogenate, yielding even more hydrogen radicals, an example of which is shown in reaction (10).Note that all these reactions are not exclusive to the Aramcomech3.0model, but are also included in the model of Wang et.al. [10] and were reported in the works of Dean [22], Roscoe and Thompson [21].Holmen [15,[23][24][25] Billaud and Guéret [14,[26][27][28][29] and Matheu et al. [30].Note that even though these reactions most likely form the dominant pathway in nonoxidative coupling of methane, their individual rates will still require tuning to accurately predict measured conversion and product selectivity in a blank system, as demonstrated in SI sections S.3 and S.10.For this reason, the data used in developing and validating the catalytic model presented in this work are provided in SI section S. 11, which can be used for tuning the individual reaction parameters in the gas-phase models.To prevent artifacts imposed by the used measurement system (i.e.reactor material, temperature profile) it can be advised to combine the data with other published experimental data [35,36,39,43].
) The catalytic reaction generating methyl-radicals, discussed at length in section 3 of the discussion and shown in SI Table S1 and SI Fig. S9, is equivalent to the combination of reactions (1) and (2) of the gas-phase mechanism discussed above.The catalytic formation of methyl radicals appears to be much faster than the slow heterogenous activation of methane in absence of catalyst according to reaction (1), in line with the observation that presence of the catalyst causes significantly increased concentration of CH 3 ⋅ radicals (Fig. 6a).Also, the concentrations of all C 2 species increase significantly in presence of catalyst (Fig. 6b and c), caused by rapid coupling of the generated methyl radicals to ethane in gas-phase (reaction (3)).The high ethylene concentration achieved inside the catalyst bed is due to the quick consecutive reactions (4) and (5).H⋅ radical concentration remains low in the catalyst bed (Fig. 6d), due to a low formation rate from C -C coupling, as well as the high reactivity of H⋅ radicals with methane, leading to a quick consumption of the produced radicals.As a reminder, chemisorbed hydrogen stays on the active site until it can desorb as a dihydrogen molecule (SI Fig. S9, [35]).Contribution of catalytic C 2 activation as well as conversion of C 2 on the reactor-wall cannot be ruled out to contribute, although our previous work [45] shows that catalytic activation of C2 is likely limited compared to gas-phase activation.
The concentration of all species presented in Fig. 6a-d for the catalyst-initiated system converge to the same concentration as observed in simulations of the empty reactor (i.e.pure gas-phase).This Fig. 7.The proposed updated autocatalytic cycle for accelerated methane activation in the gas-phase from [45], reactions numbered following the reaction scheme presented above.
effect is very rapid for all radical species (Fig. 6a, c and d) due to their high reactivity, whereas the effect is slower for ethylene (Fig. 6b).
Coke formation is omitted from all tried gas-phase models [30,44,[53][54][55][56], whereas it can significantly impact the reaction through removal of activated species from the reacting gas in the form of coke, as well as by providing a surface where activated species can adsorb on the coke, stabilizing them.Coke can either form via polymerization of hydrocarbons in the gas-phase and deposit as coke on the catalyst/reactor wall, or activated species can encounter a nucleation site on the catalyst/reactor wall leading to coke deposition [17,29].Discrimination between these two routes is difficult, but ideally the gas phase route should be included in the gas-phase models and surface deposition should be included as parameter in the surface mechanism.

Effect of C 2 addition
Ethane and ethylene addition is also known to shorten the induction phase during non-oxidative methane coupling [45,47].Experimentally, the increase in methane conversion rate was maximal at 2 % ethane or ethylene addition in a non-catalytic system [45,48].Parity plots comparing the model predictions and experimental data for an empty reactor, including the reactor temperature profile as discussed in the experimental section, have been included in SI Fig. S14a.The model accurately predicts the trends upon C 2 addition (SI Fig. S14b), although the model underpredicts the effect on methane conversion.The underprediction is caused by the stability of ethylene in the model as shown in SI Simulations of an idealized system, i.e. empty reactor completely isothermal at 1000 • C, have been performed to investigate the impact of C 2 species on the auto-catalytic cycle discussed before.Fig. 8a shows the effect of residence time on methane conversion when adding 2 vol% ethane, ethylene and acetylene.The significant effect of C 2 addition is evident when compared with the simulation for pure methane.In order to compare the effect of C 2 addition to that of catalyst initiation, in Fig. 8b the conversion curve corresponding to the simulation for pure methane is translated down by 4 % conversion and to the left by 3.3 s residence time, analogous to Fig. 5c, i.e. redefining the origin of the simulation with respect to pure-methane in Fig. 8a.The conversion axes are thus different for the pure-methane and C 2 addition curves in Fig. 8b.There is a remarkable overlap of the translated conversion curve for pure methane with the conversion curve in presence of ethane as well as ethylene.This implies that C 2 addition can effectively overcome the induction period and directly cause a reaction rate equal to pure methane in the auto-catalytic regime.Note as a qualitative observation that converting methane at 4 % corresponds stoichiometrically to the addition of 2 % C 2 hydrocarbons.
The product distributions as well as the radical concentration profiles in case of either C 2 addition or in the presence of catalyst, forming CH 3 ⋅ radicals, are identical.This makes it impossible to experimentally distinguish between catalytic CH 3 ⋅ radicals generation [35] and direct catalytic ethylene formation [43].Toraman et al. [43] show the same conclusions in simulations presented in their supporting information, where they show the concentration of the major products as function of post catalytic residence time, when either forming solely CH 3 ⋅ radicals on the catalyst or allowing for direct catalytic C 2 H 4 formation.The simulated concentrations of all products become indistinguishable after 0.5 s post catalytic residence time.This shows, together with the work presented in our publication, that the product distribution is solely determined by gas-phase reactions, the catalyst being mainly responsible for the rate of methane activation.This observation is backed-up by plotting the literature data of both empty reactors as well as catalytic reactors using the Fe©SiO 2 catalyst at higher temperatures (950+ • C) as shown in supporting information section S.10 and SI Fig. S17.Product selectivity as function of methane conversion is independent of the system used and independent of the use of the Fe©SiO 2 catalyst.This shows the need for a focus on reactor design to maximize productivity and minimize carbon formation in catalytic non-oxidative coupling of methane.The catalyst bed should, hence, be just sufficiently thick to initiate the auto-catalytic cycle, followed by a longer residence time in the hot downstream volume with minimal exposure to any surface area, e.g.reactor wall, to allow suppressing carbon formation, simultaneously maximizing methane conversion.SI section S10 and SI Fig. S17 furthermore show that more effort should be put in tuning the reaction parameters in the gas-phase models used for non-oxidative coupling of methane.It is suggested to tune either the Aramcomech3.0[53] model or the Nuigmech1.1 [54], due to their versatility in including oxygen containing reactions, or to update the model from Wang et.al. [44] from a propylene pyrolysis model to a dedicate non-oxidative methane pyrolysis model.

Conclusions
The interplay between catalytic and gas-phase reactions in nonoxidative coupling of methane was investigated.Via kinetic modelling, the considered surface reactions for a Fe©SiO 2 catalyst, involving exclusively conversion of methane to CH 3 ⋅ radicals, can accurately predict methane conversion, when used in combination with a comprehensive gas-phase reaction mechanism from literature, e.g. the AramcoMech 3.0 model.The activation barriers proposed in literature accurately describe the experimental data in the temperature range between 950 and 1100 • C. The model predicts the significant contribution of the post-catalyst volume to methane conversion.
The catalyst is proposed to enhance the formation of methyl radicals, thereby speeding up the formation of C 2 hydrocarbons via gas-phase reactions.The presence of C 2 hydrocarbons enables fast formation of The same modelled methane conversion when dosing 2 % C 2 as shown in Fig. 8a, but with the conversion curve for pure methane translated by 3.3 s residence time to the left and 4 % conversion downwards, to make the 4 % conversion point at the origin, the 2nd × and y-axis refer to the pure-methane after induction curve.All simulations done for an empty reactor system at a constant 1000 • C. H⋅ radicals, which is the key radical responsible for autocatalytic conversion of methane.The induction period observed in the gas-phase can be prevented not only by using a catalyst, but also by adding C 2 hydrocarbons in the feed, without any catalyst.All concentrations of intermediate species and products become similar to the concentrations obtained in presence or absence of catalyst, making gas-phase chemistry the sole determinant of the selectivity distribution.However, none of the tried gas-phase models accurately predicted the selectivity distribution as function of methane conversion, presented in either this work or other literature on non-oxidative coupling of methane.Further improvement of gas-phase models is required, accurately describing the product distribution in non-oxidative conversion of methane.

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Fig. 1 .
Fig. 1.(a) Three different catalyst loadings considered for the model: (1) a reactor-zone completely filled with catalyst (820 mg), from now on called 'F6 ′ fully filled for 6 cm; (2) a medium catalyst bed (260 mg) at the top of the reactor-zone, named 'T2 ′ top-filled 2 cm; (3) a medium catalyst bed (260 mg) at the bottom of the reactor zone, named 'B2 ′ bottom-filled 2 cm.(b) the differentiation the reactor modelled by the four PFRs.

Fig. 4 .
Fig. 4. (a) parity plot of model results compared to experimental results published in[41] concerning the effect of extended gas-phase residence time at 1000 • C upstream (B2) and downstream (T2) of the catalyst and (b) Modelled methane conversion for the T2 and B2 case as function of inverse total gas-flowrate; (c) &(d): parity plot comparing the selectivity of model to experimental selectivity published in[41] concerning the effect of extended gas-phase residence time at 1000 • C for the (a) B2 and (b) T2 cases; Reactor-zone at 1000 • C, pre-heater and post-heater at 400 • C; Fe©SiO 2 250-500 μm; 260 mg catalyst and 0.50 ml free volume either upstream (B2) or downstream (T2) of the catalyst.

Fig. 5 .
Fig. 5. (a) An idealized simulation of the 'T2 ′ case, uniform 1000 • C and no residence time upstream of the catalyst, achieving 4 % CH 4 conversion in the catalyst bed (b) Idealized empty reactor simulation, uniform 1000 • C (c) A comparison between the simulation shown in (a) and (b), showing a catalyst initiated system and a pure gas-phase system at a constant temperature of 1000 • C (d) shows the effect of the amount of catalyst, achieving 0 % (blank), 2 % (0.32 s catalyst residence time), 4 % (optimal, 0.48 s catalyst residence time) or 9 % (0.63 s catalyst residence time) CH 4 conversion in the catalyst bed, on the conversion curve in the post-catalytic volume.Note that the residence time on the x-axis refers only the free-volume curve the residence times for the 3 catalytic cases starts at zero at each respective left catalyst boundary; the vertical lines denote the residence in the catalyst bed.The horizontal blue dotted line denotes 4 % conversion.(For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.) Fig. S.4b as well the model not yet predicting a fast autocatalytic rate at the relatively high flowrates used for the C 2 measurements (16.6 ml⋅min − 1 , SI Fig. S3 and Fig. S10 a).Furthermore, C 2 activation on the reactor-wall may contribute to the experiments.

Fig. 8 .
Fig. 8. (a) Modelled methane conversion when dosing 2 % C 2 into the reactant mixture in a gas-phase system compared to the pure methane simulation (b)The same modelled methane conversion when dosing 2 % C 2 as shown in Fig.8a, but with the conversion curve for pure methane translated by 3.3 s residence time to the left and 4 % conversion downwards, to make the 4 % conversion point at the origin, the 2nd × and y-axis refer to the pure-methane after induction curve.All simulations done for an empty reactor system at a constant 1000 • C.