Computational investigation of acetylene hydrogenation to ethylene over transition metal–exchanged chabazite zeolites: mechanistic insights and descriptor-based predictions

Abstract

The selective hydrogenation of acetylene to ethylene is a key reaction for producing polymer-grade ethylene in industrial processes. We study the efficiency of transition metals (Fe, Co, Ni, Cu, and Zn) confined in a chabazite (CHA) zeolite for the selective hydrogenation of acetylene to ethylene using density functional theory (DFT) calculations with the M06-L functional. The proposed reaction mechanism involves H₂ dissociation, followed by the formation of an ethenyl intermediate and its subsequent conversion to ethylene. We find that Cu-CHA exhibits the highest catalytic activity among the studied catalysts and shows high selectivity for the hydrogenation of acetylene to ethylene based on the calculated overall reaction barrier and its selectivity parameter. The zeolite frameworks are found to stabilize all the species formed along the reaction pathway, but particularly at the transition states, thereby lowering all the activation barriers. We also find that the zeolite's dipole moment and metal charge are moderately accurate descriptors for predicting the overall activation barrier of the reaction. We use the sure-independence screening and sparsifying operator (SISSO) method to identify optimal nonlinear combinations of DFT-based descriptors,enabling accurate predictions without time-consuming reaction pathway calculations.

---

1. Introduction

The selective hydrogenation of acetylene to ethylene is a crucial reaction that plays an essential role in eliminating trace acetylene impurities from ethylene feedstock prior to polymerization in the polymer production process.^1–3 This reaction requires facile acetylene activation and weak ethylene adsorption to achieve high-performance catalysis. Palladium (Pd)-based catalysts exhibit high activity in acetylene hydrogenation,^4–7 but their low selectivity toward ethylene and the high cost of Pd remain challenges. Consequently, developing catalysts that combine high activity and selectivity for acetylene hydrogenation to ethylene in a cost-effective manner remains an active area of research. To this end, various catalysts have been explored, including non-noble metal catalysts^1,6,7 and metal oxide materials.^8–11

Zeolites serve as versatile silica-based supports for metal catalysts. They offer several advantages, including environmental friendliness, reusability, and regenerability, while their nanoporous frameworks impart size- and shape-selective properties. Recently, metals supported on various types of zeolites have emerged as promising catalysts for acetylene hydrogenation to ethylene such as Y-, Beta-, ZSM-12-, and CHA-supported non-noble metal catalysts.^12–24 For example, Alhashmi et al.¹⁵ reported that single-atom cobalt (Co) supported on a Y zeolite efficiently catalyzed acetylene hydrogenation, achieving an ethylene selectivity of 90% ± 2% for full acetylene conversion and exhibiting a stable performance for over 400 h and a high turnover frequency (TOF). Similarly, an iron (Fe) atom anchored on the same zeolite exhibited an ethylene selectivity of up to 93%, a high TOF and excellent stability for over 600 h.¹⁶ Chai et al.²⁰ reported the use of nickel (Ni), one of the first-row transition metals, confined within CHA zeolite frameworks, achieving 100% acetylene conversion with up to 97% selectivity toward ethylene under optimized conditions. Moreover, the introduction of modifiers further enhances the catalytic performance; for example, Sn-modified Ni/ZSM-12¹⁹ and alkali metal cation–modified Ni–CHA²⁰ achieved complete acetylene conversion with high ethylene selectivity. These findings motivated us to investigate and compare the activity and selectivity of non-noble transition metals supported on zeolites, along with the influence of zeolite frameworks on acetylene hydrogenation to ethylene, which has not yet been reported.

To gain a fundamental understanding of the behaviors of zeolites, computational methods, especially those based on density functional theory (DFT) calculations, provide atomistic insights into energetics, electronic structures, and molecular geometries across the potential energy surface, which is information often inaccessible by experimental techniques. They especially enable the identification of active sites, adsorption characteristics and reaction pathways, offering valuable guidance for the rational design of efficient and selective catalysts. Among these quantities, the activation barrier is particularly crucial for elucidating trends in the catalytic performance of zeolite-based systems. Its determination requires locating transition-state structures, which can be computationally demanding.²⁵ Therefore, more cost-effective approaches, particularly those identifying DFT-based descriptors, are valuable for estimating activation barriers. Combinations of descriptors can achieve good predictive accuracy by using the sure independence screening and sparsifying operator (SISSO) method²⁶ if overfitting is carefully avoided. Recently, this approach has been successfully applied not only to predict activation energies but also to estimate transition-state energies for reactions on porous materials, including zeolites and metal–organic frameworks.^27–30 To the best of our knowledge, the systematic application of SISSO to acetylene hydrogenation to ethylene over metal sites confined within zeolite frameworks has not previously been reported.

In this work, we investigate the catalytic activity and selectivity of transition metals (Fe, Co, Ni, Cu, and Zn) confined within the CHA zeolite framework for acetylene hydrogenation to ethylene using DFT-based quantum chemical calculations, providing detailed mechanistic and energetic insights. The effects of different metal species on catalytic performance and ethylene selectivity are systematically evaluated. Additionally, several DFT-based descriptors related to the geometric and electronic properties of the zeolites are identified to predict the overall activation barrier. The SISSO method is employed to derive compact nonlinear expressions.

2. Models and methods

The chabazite (CHA) zeolite is represented by a cluster model consisting of 42 tetrahedral (T) atoms, encompassing both double six-membered rings (D6R) and the CHA cage, as illustrated in Fig. 1. This model is derived from the crystal lattice of the CHA zeolite. The Si/Al ratio is 20, consistent with previously reported computational investigations of Ni-exchanged CHA catalysts.²⁰ The metal-exchanged zeolite is constructed by substituting two silicon atoms in a six-membered ring by aluminum atoms and a transition metal atom (Fe, Co, Ni, Cu, or Zn) is introduced at the six-membered ring site.^20,27 In this configuration, the metal center is four-coordinated by the framework oxygen atoms, which is consistent with X-ray absorption spectroscopy (XAS) measurements reporting a coordination number of approximately 4.0 for the Ni-CHA case.^20,31

Fig. 1 Cluster model of the metal-supported CHA zeolite used in this work. Si, Al, O, C, H and metal atoms are represented in green, pink, red, gray, white, and blue, respectively.

All geometry optimizations of reactants, intermediates, products, and transition states in this work were carried out using the M06-L density functional³² as implemented in Gaussian 16.³³ This functional incorporates medium-range dispersion interactions within its parameterization and is well-suited for systems containing transition metals.³⁰ Its robustness and accuracy in describing adsorption and reaction mechanisms on zeolite frameworks have been established in our previous studies.^34–41 The def2-SVP basis set was used to describe the C, O, and H atoms, while the Stuttgart–Dresden pseudopotential with a double-ζ basis set was applied for the transition metals. During geometry optimization, both the interacting molecules and the zeolite framework were fully relaxed with only the hydrogen atoms used for cluster termination being fixed. To characterize each stationary point along the reaction pathway, vibrational frequency analyses were conducted on the optimized structures at the same level of theory. These calculations also provide thermodynamic properties. They were evaluated at 453.15 K, corresponding to experimental studies on acetylene hydrogenation over zeolites.²⁰ The Gibbs free energy (G) is evaluated according to eqn (1), as follows:

G = E_el + G_corr (1)

where E_el symbolizes the total electronic energy and G_corr denotes the free energy correction. G_corr is gained from eqn (2), as follows:

G_corr = E_tot + k_BT − TS_tot (2)

where E_tot denotes the energy contribution from the internal degrees of freedom, k_B is the Boltzmann constant, T is temperature and S_tot is the internal entropy. The values for E_tot and S_tot contain translational, rotational, and vibrational energy contributions. The relative free energy along the reaction coordinates is defined as follows:

ΔG = G_complex − G_probe – G_CHA (3)

where G_complex represents the total Gibbs free energy of the adsorption state, transition state, intermediate, and product; G_probe is the Gibbs free energy of the isolated reactants, and G_CHA is the Gibbs free energy of the isolated CHA zeolite.

3. Results and discussion

3.1. Reaction mechanisms and energetic trends

Fig. 1 shows the cluster model of M-CHA. As mentioned above, the transition-metal atom (Fe, Co, Ni, Cu, or Zn) is located at the six-membered ring site of the CHA framework. In this configuration, the metal center is fourfold coordinated by framework oxygen atoms. This coordination environment is in good agreement with previous experimental characterizations and computational studies, particularly for Ni- and Cu-exchanged CHA systems, which report metal coordination numbers of approximately four.^20,42–44 Moreover, in the Cu case, the calculated nearest-neighbor Cu–O bond distance is 1.93 Å, which is in good agreement with our previously reported data for Cu-CHA (1.95 Å).⁴² The consistency with these earlier studies supports the reliability of the present cluster model in capturing the local geometric and electronic environment of the active metal centers.

We further discuss the reaction mechanism and the energetic trends of acetylene hydrogenation to ethylene on M-CHA (M = Zn, Cu, Ni, Co, or Fe). The reaction proceeds in three main steps: (i) hydrogen molecule dissociation, (ii) the formation of a surface ethenyl intermediate, and (iii) ethylene formation, as shown in Fig. 2.²⁰ The relative free energy profile is given in Table 1, and the corresponding Gibbs free energy barriers (ΔG^#) are illustrated in Fig. 3. The optimized structures for the reaction on M-CHA are displayed in Fig. S1–S5 of the SI.

Fig. 2 Reaction mechanism for the acetylene hydrogenation to ethylene on M-CHA. Si, Al, O, C, H and metal atoms are represented in green, pink, red, gray, white, and blue, respectively.

Table 1 Relative Gibbs free energies for the hydrogenation of acetylene to ethylene on M-CHA

Species	Relative Gibbs free energy (kcal mol^−1)
	Fe-CHA	Co-CHA	Ni-CHA	Cu-CHA	Zn-CHA
Ads	10.4	3.2	3.3	6.7	5.1
TS1	32.7	24.6	17.7	18.7	21.0
Int1	33.1	25.0	18.1	16.6	13.2
Int1-Ace	30.3	23.4	14.8	13.4	16.9
TS2	40.1	35.9	25.0	22.7	43.0
Int2	−3.0	−11.0	−16.7	−20.1	−19.1
TS3	−1.5	−7.8	−13.5	−14.6	−13.5
Prod	−42.4	−49.1	−49.6	−39.6	−46.3

---

Fig. 3 Gibbs free energy barrier for the acetylene hydrogenation to ethylene on M-CHA calculated at 453.15 K.

The reaction starts with H₂ binding to the metal site (Ads) in a symmetric configuration, with identical H⋯M distances. Upon adsorption, the H₂ molecule becomes more positive, while the positive charge on the metal site decreases relative to the isolated zeolite systems, indicating electron transfer from H₂ to the zeolite active site. This charge transfer leads to elongation of the H–H bond from 0.75 to 0.76–0.78 Å. In line with this minor geometrical change, the adsorption free energy of H₂ adsorption is endothermic at 453.15 K for all the systems (see Table 1). The adsorbed H₂ then dissociates via the transition state (TS1), in which the H–H bond is cleaved and the two hydrogen atoms are simultaneously transferred to the metal site and the oxygen atom (O1) of the zeolite, as illustrated in Fig. 2. In this state, the H–H bond distance is elongated, while the M⋯H and O1⋯H distances are shortened, indicating bond formation. The dissociation of H₂ proceeds via heterolytic cleavage. This is shown by the NBO partial charge analysis, where H1 gains a negative partial charge while H2 acquires a positive one. The activation barriers for this step (Ea1), relative to the energy of the isolated molecule, increase as follows: Cu-CHA (12.0 kcal mol⁻¹) < Ni-CHA (14.4 kcal mol⁻¹) < Zn-CHA (15.9 kcal mol⁻¹) < Co-CHA (21.4 kcal mol⁻¹) < Fe-CHA (22.3 kcal mol⁻¹), as shown in Fig. 3. These activation energies are in the same range as those previously reported⁴⁵ for H–H activation on metal alkoxide MOFs. The value obtained for Ni-CHA is also consistent with previous theoretical studies.²⁰ The lower barrier observed for Cu-CHA can be ascribed to its shorter H–H bond distance in the TS1 structure, which facilitates earlier bond activation, and thus stabilizes the transition state. The NBO charge analysis also shows that the hydride species in TS1 of Cu-CHA is less negative than the ones associated with the other metals (see Table S1). This stems from a larger donation of electron density from the H₂ σ-bonding orbital into the vacant orbitals of Cu, which contributes to the lower activation barrier. The transition state subsequently leads to the formation of a metal-bound hydride and a proton bound to the zeolite oxygen atom (Int1).

An acetylene molecule diffuses and adsorbs on the metal-bound hydride site (Int1-Ace) with calculated relative energies, referenced to the isolated reactants, ranging from 13.4 to 30.3 kcal mol⁻¹ (see Table 1). Acetylene is then hydrogenated to form the ethenyl intermediate via the TS2 transition state. The M–H1 bond breaks, transferring H1 to the C1 carbon of acetylene, while a new bond forms between C2 and metal. TS2 is considered a four-center transition state, formed through the breaking of the M–H1 bond and the formation of the C1–H1 and M–C2 bonds. The activation barriers (Ea2) relative to the co-adsorption complex are 9.8, 12.5, 10.2, 9.3, and 26.1 kcal mol⁻¹ for Fe-, Co-, Ni-, Cu-, and Zn-CHA, respectively (see Fig. 3). Next, the ethenyl intermediate (Int2) is formed at the active site of M-CHA. This species has been experimentally observed at 453 K on the Ni-CHA zeolite.²⁰ As revealed by the NBO charge analysis (see Table S1), the lower activation barrier for Cu-CHA originates from the smaller positive charge on the Cu center in Int1-Ace compared to the other metals. This causes a weaker electrostatic interaction with the hydride species and a lower activation barrier for M–H bond cleavage. The complexation energies for the formation of this intermediate are exothermic in all the systems (−3.0 to −20.1 kcal mol⁻¹), as shown in Fig. 3.

Following the formation of the ethenyl intermediate, the reaction proceeds through its hydrogenation, resulting in the formation of ethylene. The transition state of this step, as shown in Fig. 2 (TS3), involves the concerted cleavage of the O1–H2 bond and formation of the H2–C2 bond, resulting in an increased O1–H2 distance and a decreased H2⋯C2 distance. The activation barriers (Ea3) in this step are low, ranging from 1.5 kcal mol⁻¹ for Fe-CHA to 5.6 kcal mol⁻¹ for Zn-CHA. This indicates that ethenyl can be easily converted to the ethylene product. The lower activation barrier for Fe-CHA can be attributed to the less negative charge on the Fe center in TS3 (see Table S1). Thus, proton release is facilitated and the energy barrier is lowered. Subsequently, the ethylene product (Prod) adsorbs on the Ni active site through its π bond, with a complexation energy ranging from −39.6 to −49.6 kcal mol⁻¹ relative to the isolated molecule. Finally, the desorption of the product from the M-CHA surface requires 3.9, 10.6, 11.1, 1.1, and 7.8 kcal mol⁻¹ for Fe-, Co-, Ni-, Cu-, and Zn-CHA, respectively.

3.2. Effect of transition metals

As shown in Fig. 3, for Fe-, Co-, Ni-, and Cu-CHA, the activation free energy for H–H bond dissociation is higher than that for either of the hydrogenation steps for acetylene. This is in contrast with Zn-CHA, where the first hydrogenation of acetylene requires the highest barrier among the steps. However, for all systems, the H₂ adsorption complex (Ads) is identified as the rate-determining intermediate, with transition state TS2 leading to the formation of the ethenyl intermediate as the rate-determining transition state. Therefore, the overall free energy barriers (ΔG^#_overall) are calculated from the energy differences between these two states, as defined below:

ΔG^#_overall = G_TS2 − G_Ads

The ΔG^#_overall values are 29.7, 32.7, 21.8, 16.0, and 37.9 kcal mol⁻¹ for Fe-, Co-, Ni-, Cu-, and Zn-CHA, respectively (see Fig. 4). The lower overall activation barrier for Cu-CHA is associated with the higher energetic stability of TS2 (see Table 1). This enhanced stability of TS2 for Cu-CHA is directly linked to the negligible charge accumulation on the ethenyl intermediate, which remains close to neutral (see Table S1). The minimized internal electrostatic repulsion results in the lowest transition-state energy among the studied metal centers.

Fig. 4 Overall Gibbs free energy barrier (ΔG^#_overall) and selectivity parameter (ΔΔG) for the acetylene hydrogenation to ethylene on M-CHA calculated at 453.15 K.

To check the importance of dispersion effects, we compared the M06-L results with those obtained using M06-L-D3, a model that explicitly includes dispersion terms, for ΔG^#_overall. The results show that the trends in ΔG^#_overall across the metal series remain unchanged (Table S2). This indicates that the inclusion of dispersion corrections is not crucial for our conclusions.

We further validated the computational methodology employed in this work by benchmarking the M06-L functional against another commonly used functional, B3LYP combined with both Grimme's D3 dispersion correction and the D3(BJ) damping scheme, and by comparing results obtained with the def2-SVP and def2-TZVP basis sets (see Tables S2 and S3), respectively. In all cases, the ΔG^#_overall values exhibit the same trends as found with M06-L, identifying Cu-CHA as the most favorable candidate. Therefore, we believe that M06-L in combination with the chosen level of theory provides a reliable description of the systems studied.

In addition to the catalytic activity, the selectivity for ethylene formation represents a crucial factor in acetylene hydrogenation. As reported previously, ethylene desorption from the catalyst surface is generally favored over its further hydrogenation.^46,47 The latter requires the dissociation of an additional H₂ molecule as the initial step for hydrogenating the ethylene formed in the primary reaction. Consequently, high selectivity towards ethylene formation is achieved when the free energy barrier for ethylene desorption is lower than that for a subsequent H₂ dissociation step, leading to ethylene over-hydrogenation. The difference between these two barriers has been used as selectivity parameter in previous studies on selective acetylene hydrogenation on metal–alkoxide MOFs.³⁸ Accordingly, eqn (4) defines this difference between two relative free energies:

ΔΔG = ΔG_over – ΔG_des (4)

where ΔG_over is the activation free energy of H₂ dissociation in the presence of adsorbed ethylene⁴⁵ and ΔG_des is the free energy barrier for ethylene desorption. A larger positive ΔΔG value corresponds to higher selectivity toward ethylene formation. The optimized structures of H₂ dissociation in the presence of adsorbed ethylene are given in Fig. S6 of the Supporting Information. As shown in Fig. 4, all the systems exhibit positive ΔΔG values at the studied temperature, indicating that ethylene desorption (ΔG_des) is energetically more favorable than further hydrogenation (ΔG_over). This suggests that all the catalysts possess high selectivity toward ethylene formation, with Fe-CHA exhibiting the highest selectivity, as reflected by its more positive ΔΔG value compared to the others. In the case of Ni-CHA, this result is consistent with previous experimental studies reporting its high selectivity for acetylene hydrogenation to ethylene.²⁰ Combining ΔΔG with the ΔG^#_overall values reported above, Cu-CHA exhibits the highest catalytic activity among the studied catalysts, together with high selectivity for acetylene hydrogenation to ethylene. It remains the most promising candidate as a catalyst for this reaction, also considering its selectivity.

Concerning the issue of ethylene selectivity discussed above, we explicitly investigated the possible over-hydrogenation of ethylene to ethane on all the studied zeolites. This process proceeds via two consecutive hydrogenation steps: firstly, hydrogenation of ethylene to form the ethyl intermediate, followed by a second hydrogenation step leading to ethane formation (see Fig. S7). The relative free energies are reported in Table S4. The calculated overall free-energy barriers for the hydrogenation of ethylene to ethane, defined as the free-energy difference between the co-adsorbed ethylene and the H₂ state (Int3; rate-determining intermediate) and the transition state for ethyl formation (TS5; rate-determining transition state) are 45.1, 43.3, 31.0, 16.1, and 50.7 kcal mol⁻¹ for Fe-, Co-, Ni-, Cu-, and Zn-CHA, respectively. These barriers are higher than those associated with the initial acetylene hydrogenation to ethylene and also exceed the free energy barriers for ethylene desorption reported above. These results indicate that once formed, ethylene preferentially desorbs from the catalyst surface rather than undergoing further hydrogenation to ethane. This kinetic preference supports the notion of high selectivity toward ethylene formation in these systems.

To reveal the influence of the zeolite framework, we carried out single-point calculations on the 42T geometries using a reduced 12T cluster model of M-CHA, which preserves only the D6R structure (see Fig. S8). The relative energies of the systems along the reaction coordinate are summarized in Table S1. We observed a strong linear correlation between the relative energies calculated from the 42T and 12T clusters, as indicated by the coefficients of determination (R²) ranging from 0.928 to 0.979, demonstrating that both models follow the same trend. Nevertheless, the relative energies obtained from the 42T cluster model are consistently lower than those from the 12T cluster, particularly for the reaction intermediates and transition states. The energy differences fall in the ranges of 6.3–15.7 kcal mol⁻¹ for Fe-CHA, 0.9–11.0 kcal mol⁻¹ for Co-CHA, 2.4–13.6 kcal mol⁻¹ for Ni-CHA, 0.3–15.9 kcal mol⁻¹ for Cu-CHA and 4.6–14.7 kcal mol⁻¹ for Zn-CHA, as summarized in Table S5. These results reflect the stabilizing effect of the extended zeolite framework. Fig. 5 shows the result of the linear fit between the ΔG^#_overall from the large and the small cluster models. The R² value of this regression is 0.929. The trends in ΔG^#_overall are the same for all the metal systems when the 12T model is used. However, the ΔG^#_overall values for the 12T cluster are consistently higher by 11.5–14.3 kcal mol⁻¹ than those from the 42T cluster for all the zeolites (see Table S5). This indicates that the extended zeolite framework stabilizes the species formed along the reaction coordinates and reduces the overall barrier, thereby enhancing the activity of the metal active site.

Fig. 5 Linear scaling relationship between the overall activation barriers obtained from the 12T and 42T models.

To ensure that the 42T cluster model used in this work is converged with respect to cluster size, single-point calculations of ΔG^#_overall were performed at the same level of theory on extended 78T cluster models based on the structures optimized with the 42T model (see Fig. S9). The ΔG^#_overall values of all the systems using the 42T model are nearly identical to those obtained with the 78T model, and more importantly, the trends of ΔG^#_overall across the metal series remain unchanged (see Table S6). These results indicate that the practical 42T cluster model is sufficiently large to accurately represent the interactions between the adsorbates and the zeolite framework.

3.3. Descriptors

We attempt to develop descriptors capable of estimating the overall activation barrier without the need for explicit TS searches, which are computationally very demanding and rely on prior chemical insight. For this purpose, reaction data obtained with and without the zeolite framework are utilized. We first correlate the overall activation barriers (ΔG^#_overall) with single descriptors derived from the geometric properties of the zeolites, including the average M–O bond distance between the metal center and its three nearest oxygen atoms (B_M–O; M = Fe, Co, Ni, Cu and Zn), dipole moment (Di), highest occupied molecular orbital (HOMO; E_H) energy, lowest unoccupied molecular orbital (LUMO; E_L) energy, electronic chemical potentials (μ), global electronic hardness indices (η), global electrophilicity indices (ω), and Mulliken charge of the metal center (ρ).

The R² and root-mean-square errors (RMSE) for the linear correlations between the descriptors and ΔG^#_overall are presented in Table S7. A reasonably strong linear relationship is observed for Di and ρ, with R² values of 0.718 and 0.817 and corresponding RMSEs of 5.3 and 4.3 kcal mol⁻¹, respectively. ΔG^#_overall decreases with decreasing metal charge and increasing dipole moment. These linear relationships for predicting ΔG^#_overall are:

ΔG^#_overall = −12.065 × Di + 69.712 kcal mol⁻¹ (5)

ΔG^#_overall = 72.799 × ρ + 2.213 kcal mol⁻¹ (6)

The predictive power of these simple regression models remains moderate, as reflected by the R² and RMSE values reported above. To enhance it, we employ the SISSO method and generate more complex descriptors by screening a large set of candidates formed from combinations of primary descriptors through various mathematical expressions. In this work, 8 primary descriptors, as outlined above, are used as inputs. These descriptors are expanded by applying various mathematical operators (+, −, ×, ÷, exp, exp−, ⁻¹, ², ³, √, and | |) to generate nonlinear expressions forming the candidate descriptor space. Due to the limited size of the dataset, all 10 data points were used for model training, which is a common and accepted practice in SISSO applications involving small datasets. Within the SISSO framework, the Sure Independence Screening (SIS) procedure is first applied to rank candidate descriptors according to their correlation with the target property, ΔG^#_overall. Subsequently, the sparsifying operator based on L₀ regularization is used to identify the most relevant low-dimensional descriptors. The descriptor dimensionality (D) is limited to 3, with a maximum of 3 operations per term, resulting in a total of 3 model levels.

Table 2 lists the comprehensive mathematical equations derived for each model complexity. Fig. 6 shows a comparison between the ΔG^#_overall values from DFT and those predicted by SISSO, including the corresponding R² and RMSE for the trained models. The increase in the coefficient of determination, R², from 0.938 for the 1D model to 0.975 and 0.997 for the 2D and 3D models, respectively, shows that the reproduction of the energetic trends improves progressively. The corresponding RMSE values decrease from 2.5 to 1.6 and 0.5 kcal mol⁻¹, reflecting a substantial reduction in the errors in predicting ΔG^#_overall for the higher-dimensional models. The robustness of these predictive models is validated to prevent overfitting using the leave-one-group-out cross validation (LOOCV) procedure.⁴⁸ The R² and RMSE values are used to evaluate the overall cross-validation performance. The high R² values obtained (0.938–0.997) together with low RMSEs (0.1–1.0 kcal mol⁻¹) indicate that the SISSO models for predicting ΔG^#_overall are robust and not prone to overfitting. It should be noted that the SISSO analysis in this work is based on a limited but chemically consistent dataset comprising five transition-metal-exchanged CHA zeolites described by two cluster models. SISSO has been shown to be effective in small-data regimes for descriptor discovery.⁴⁹ Nevertheless, expanding the dataset to include additional metals and different zeolite topologies would be beneficial for assessing both predictive power and descriptor transferability. These extensions will be the focus of future work.

Table 2 Best SISSO models for predicting the overall activation barriers (ΔG^#_overall)

Descriptor dimension	SISSO model
1D
2D
3D

---

Fig. 6 Comparison of ΔG^#_overall values obtained from DFT and SISSO in one-dimensional (1D) (a), two-dimensional (2D) (b), and three-dimensional (3D) (c) models, and the corresponding R² and RMSE values (d).

We can try to rationalize the physical meaning of the nonlinear descriptors selected by the SISSO algorithm for describing the overall activation free energy. The 1D descriptor captures the dominant contribution to the activation barrier by combining the average metal–oxygen bond distance (B_M–O) with the dipole moment (Di) and the HOMO energy (E_H). The inverse cubic dependence on B_M–O highlights the dominant role of the metal–oxygen interaction strength in governing the overall activation barrier, whereby shorter B_M–O distances tend to correspond to lower barriers. This geometric contribution is further modulated by the Di/E_H term, which reflects the influence of electronic polarization and donor ability of the metal-exchanged CHA zeolites. In particular, a higher dipole moment is associated with a reduced activation barrier. Overall, the 1D descriptor containing three basis properties captures the leading geometric–electronic control of the reaction barrier and provides a physically transparent baseline model for describing the catalytic activity trend.

The 2D descriptor extends the 1D model by adding an additive term that contains the Mulliken charge of the metal center (ρ), coupled with the mismatch between B_M–O and Di effects. This additional term accounts for local charge redistribution at the active site, indicating that the activation barrier is not determined solely by the bond strength but is also sensitive to the degree of charge accumulation and polarization around the metal center. The improved performance of the 2D model demonstrates the importance of simultaneously considering geometric structure and electronic charge effects.

The 3D descriptor further includes higher-order electronic properties, including the LUMO energy (E_L), global hardness (η), and electrophilicity index (ω), together with charge ρ. The inclusion of these electronic terms further improves the prediction accuracy by approximately the same amount as going from the one- to the two-dimensional expression.

4. Conclusions

The hydrogenation of acetylene to ethylene over different metal atoms (Fe, Co, Ni, Cu and Zn) incorporated in a CHA zeolite framework has been investigated by DFT calculations. A three-step mechanism for the hydrogenation of acetylene to ethylene is proposed and investigated, including (i) H₂ dissociation, (ii) the formation of an ethenyl intermediate, and (iii) ethylene product formation. The rate-determining intermediate is the adsorbed H₂ species (Ads), while the rate-determining transition state (TS2) corresponds to the formation of the ethenyl intermediate. The overall activation free energies are calculated to be 29.7, 32.7, 21.8, 16.0, and 37.9 kcal mol⁻¹ for Fe-, Co-, Ni-, Cu-, and Zn-CHA, respectively. Ethylene desorption is energetically favorable against H₂ dissociation in the presence of adsorbed ethylene, as indicated by the lower desorption barrier. This suggests that all the zeolite systems enhance the selectivity toward the desired ethylene product. Based on the results, Cu-CHA has the highest catalytic activity. Its selectivity is also high, highlighting its potential as a catalyst for this reaction. We show that the extended zeolite framework not only plays a crucial role in stabilizing all the species along the reaction pathway, but also in reducing the overall activation barrier. Based on the linear regression analysis using single descriptors, the dipole moment of the zeolite and the partial charge of the metal site are moderately effective descriptors for predicting the overall activation barriers. The predictive accuracy can be enhanced using the SISSO statistical learning approach, which systematically combines DFT-derived descriptors from isolated zeolites to identify the most relevant predictors. The resulting compact expressions have been checked against overfitting and enable an accurate estimation of the overall activation barriers, allowing to bypass the time-consuming reaction pathway calculations in catalyst screening in many cases.

Author contributions

Thana Maihom: investigation, conceptualization, data curation, writing – original draft, writing – review and editing. Michael Probst: writing – review and editing. Jumras Limtrakul: supervision.

Conflicts of interest

There are no conflicts to declare.

Data availability

The Cartesian coordinates of all the optimized structures can be found at https://doi.org/10.5281/zenodo.18419313.

Supplementary information (SI) is available. See DOI: https://doi.org/10.1039/d5cp04172d.

Acknowledgements

This work was supported by the Kasetsart University Research and Development Institute, KURDI (FF(KUSRIU)5.67).

Notes and references

1. F. Studt, F. Abild-Pedersen, T. Bligaard, R. Z. Sørensen, C. H. Christensen and J. K. Nørskov, Science, 2008, 320(5881), 1320–1322.
2. M. T. Ravanchi, S. Sahebdelfar and S. Komeili, Rev. Chem. Eng, 2018, 34(2), 215–237.
3. L. Zhang, M. Zhou, A. Wang and T. Zhang, Chem. Rev., 2020, 120(2), 683–733.
4. C. W. Chan, A. H. Mahadi, M. M. Li, E. C. Corbos, C. Tang, G. Jones, W. C. Kuo, J. Cookson, C. M. Brown, P. T. Bishop and S. C. Tsang, Nat. Commun., 2014, 5, 5787.
5. L. B. Ding, H. Yi, W. H. Zhang, R. You, T. Cao, J. L. Yang, J. L. Lu and W. X. Huang, ACS Catal., 2016, 6, 3700–3707.
6. M. Armbrüster, K. Kovnir, M. Friedrich, D. Teschner, G. Wowsnick, M. Hahne, P. Gille, L. Szentmiklósi, M. Feuerbacher, M. Heggen, F. Girgsdies, D. Rosenthal, R. Schlögl and Y. Grin, Nat. Mater., 2012, 11, 690–693.
7. Y. Liu, X. Liu, Q. Feng, D. He, L. Zhang, C. Lian, R. Shen, G. Zhao, Y. Ji, D. Wang, G. Zhou and Y. Li, Adv. Mater., 2016, 28, 4747–4754.
8. G. Vilé, P. Dähler, J. Vecchietti, M. Baltanás, S. Collins, M. Calatayud, A. Bonivardi and J. Pérez-Ramírez, J. Catal., 2015, 324, 69–78.
9. K. Werner, X. Weng, F. Calaza, M. Sterrer, T. Kropp, J. Paier, J. Sauer, M. Wilde, K. Fukutani, S. Shaikhutdinov and H. J. Freund, J. Am. Chem. Soc., 2017, 139, 17608–17616.
10. M. Tejeda-Serrano, J. R. Cabrero-Antonino, V. MainarRuiz, M. Lopez-Haro, J. C. Hernandez-Garrido, J. J. Calvino, A. Leyva-Perez and A. Corma, ACS Catal., 2017, 7, 3721–3729.
11. S. Zhou, L. Gao, F. Wei, S. Lin and H. Guo, J. Catal., 2019, 375, 410–418.
12. J. Vito and M. Shetty, ACS Appl. Mater. Interfaces, 2024, 16, 67010–67027.
13. X. Deng, R. Bai, Y. Chai, Z. Hu, N. Guan and L. Li, CCS Chem., 2022, 4(3), 949–962.
14. Y. X. Zhao, Ö. Bozkurt, S. F. Kurtoglu-Öztulum, M. S. Yordanli, A. S. Hoffman, J. Y. Hong, J. E. Perez-Aguilar, A. Saltuk, D. Akgül, O. Demircan, T. A. Atesin, V. Aviyente, B. C. Gates, S. R. Bare and A. Uzun, J. Catal., 2024, 429, 115196.
15. E. Alhashmi, G. Ebri and K. Hellgardt, ACS Catal., 2025, 15, 4121–4140.
16. E. Alhashmi, G. Ebri and K. Hellgardt, Chem. Eng. J., 2025, 524, 168933.
17. W. Huang, J. R. McCormick, R. F. Lobo and J. G. Chen, J. Catal., 2007, 246(1), 40–51.
18. S. Hu, C. Zhang, M. Wu, R. Ye, D. Shi, M. Li, P. Zhao, R. Zhang and G. Feng, Catalysts, 2022, 12(9), 1072.
19. D. Wang, R. Ye, C. Zhang, C. Jin, Z. Lu, M. Shakouri, B. Han, T. Wang, Y. Zhang, R. Zhang, Y. Hu, J. Zhou and G. Feng, Energy Fuels, 2023, 37, 13305–13318.
20. Y. Chai, G. Wu, X. Liu, Y. Ren, W. Dai, C. Wang, Z. Xie, N. Guan and L. Li, J. Am. Chem. Soc., 2019, 141, 9920–9927.
21. J. Vito and M. Shetty, ACS Appl. Mater. Interfaces, 2024, 16(49), 67010–67027.
22. Y. Zhao, Ö. D. Bozkurt, S. F. Kurtoğlu-Öztulum, M. S. Yordanli, A. S. Hoffman, J. Hong, J. E. Perez-Aguilar, A. Saltuk, D. Akgül, O. Demircan, T. A. Ateşin, V. Aviyente, B. C. Gates, S. R. Bare and A. Uzun, J. Catal., 2024, 429, 115196.
23. M. S. Yordanlı, R. Escobar, J. Meza, D. Akgül, Y. Zhao, A. Uzun, F. Ahu Akin, V. Aviyente, A. C. Atesin and T. A. Ateşin, ChemPhysChem, 2025, 26, e202400867.
24. E. Alhashmi, G. Ebri and K. Hellgardt, Adv. Sustainable Syst., 2025, 9, e01017.
25. P. Y. Ayala and H. B. Schlegel, J. Chem. Phys., 1997, 107(2), 375–384.
26. R. Ouyang, S. Curtarolo, E. Ahmetcik, M. Scheffler and L. M. Ghiringhelli, Phys. Rev. Mater., 2018, 2, 083802.
27. H. A. Doan, X. Wang and R. Q. Snurr, J. Phys. Chem. Lett., 2023, 14, 5018–5024.
28. T. Maihom, J. Sittiwong, M. Probst, B. Boekfa, C. Wattanakit and J. Limtrakul, J. Catal., 2024, 437, 115656.
29. X. Wang, K. Shi, A. Peng and R. Q. Snurr, ACS Catal., 2024, 14(24), 18708–18721.
30. T. Maihom, B. Boekfa, S. Boonya-Udtayan, P. Ketwong, C. Wattanakit, M. Probst and J. Limtrakul, J. Phys. Chem. C, 2025, 129(10), 5182–5190.
31. X. Deng, B. Qin, R. Liu, X. Qin, W. Dai, G. Wu, N. Guan, D. Ma and L. Li, J. Am. Chem. Soc., 2021, 143, 20898–20906.
32. Y. Zhao and D. G. Truhlar, Acc. Chem. Res., 2008, 41, 157–167.
33. M. J. Frisch, G. W. Trucks, H. B. Schlegel, G. E. Scuseria, M. A. Robb, J. R. Cheeseman, G. Scalmani, V. Barone, G. A. Petersson, H. Nakatsuji, X. Li, M. Caricato, A. V. Marenich, J. Bloino, B. G. Janesko, R. Gomperts, B. Mennucci, H. P. Hratchian, J. V. Ortiz, A. F. Izmaylov, J. L. Sonnenberg, D. Williams-Young, F. Ding, F. Lipparini, F. Egidi, J. Goings, B. Peng, A. Petrone, T. Henderson, D. Ranasinghe, V. G. Zakrzewski, J. Gao, N. Rega, G. Zheng, W. Liang, M. Hada, M. Ehara, K. Toyota, R. Fukuda, J. Hasegawa, M. Ishida, T. Nakajima, Y. Honda, O. Kitao, H. Nakai, T. Vreven, K. Throssell, J. A. Montgomery Jr., J. E. Peralta, F. Ogliaro, M. J. Bearpark, J. J. Heyd, E. N. Brothers, K. N. Kudin, V. N. Staroverov, T. A. Keith, R. Kobayashi, J. Normand, K. Raghavachari, A. P. Rendell, J. C. Burant, S. S. Iyengar, J. Tomasi, M. Cossi, J. M. Millam, M. Klene, C. Adamo, R. Cammi, J. W. Ochterski, R. L. Martin, K. Morokuma, O. Farkas, J. B. Foresman and D. J. Fox, Gaussian 16, Revision B.01, Gaussian, Inc., Wallingford CT, 2016.
34. Y. Zhao and D. G. Truhlar, J. Phys. Chem. C, 2008, 112, 6860–6868.
35. T. Maihom, B. Boekfa, J. Sirijaraensre, T. Nanok, M. Probst and J. Limtrakul, J. Phys. Chem. C, 2009, 113, 6654–6662.
36. T. Maihom, P. Pantu, C. Tachakritikul, M. Probst and J. Limtrakul, J. Phys. Chem. C, 2010, 114, 7850–7856.
37. T. Maihom, M. Probst and J. Limtrakul, J. Phys. Chem. C, 2014, 118, 18564–18572.
38. Y. Injongkol, T. Maihom, P. Treesukul, J. Sirijaraensre, B. Boekfa and J. Limtrakul, Phys. Chem. Chem. Phys., 2017, 19, 24042–24048.
39. A. Prasertsab, T. Maihom, M. Probst, C. Wattanakit and J. Limtrakul, Inorg. Chem., 2018, 57, 6599–6605.
40. W. Thongnuam, T. Maihom, S. Choomwattana, Y. Injongkol, B. Boekfa, P. Treesukol and J. Limtrakul, Phys. Chem. Chem. Phys., 2018, 20, 25179–25185.
41. W. Jeevapong, J. Sittiwong, M. Probst, B. Boekfa, C. Wattanakit, T. Maihom and J. Limtrakul, J. Phys. Chem. C, 2023, 127, 8473–8481.
42. U. Deka, A. Juhin, E. A. Eilertsen, H. Emerich, M. A. Green, S. T. Korhonen, B. M. Weckhuysen and A. M. Beale, J. Phys. Chem. C, 2012, 116, 4809–4818.
43. A. Godiksen, F. N. Stappen, P. N. R. Vennestrøm, F. Giordanino, S. B. Rasmussen, L. F. Lundegaard and S. Mossin, J. Phys. Chem. C, 2014, 118, 23126–23138.
44. M. Moreno-Gonzalez, B. Hueso, M. Boronat, T. Blasco and A. Corma, Ammonia-containing species formed in Cu-Chabazite as per in situ EPR, solid-state NMR, and DFT calculations, J. Phys. Chem. Lett., 2015, 6, 1011–1017.
45. A. Ganai and P. Sarkar, J. Phys. Chem. C, 2024, 128, 7913–7925.
46. H. Feng, H. Ding, S. Wang, Y. Liang, Y. Deng, Y. Yang, M. Wei and X. Zhang, ACS Appl. Mater. Interfaces, 2022, 14, 25288–25296.
47. Y. Jie, J. Y. Guo, H. Lu, T. Y. Shen, Y. F. Zhang, J. Zhao, C. X. Hu, M. Pu and H. Yan, J. Phys. Chem. C, 2025, 129(29), 13219–13231.
48. Y. Jung and J. A. Hu, J. Nonparametric Stat., 2015, 27, 167–179.
49. M. Khatamirad, E. Fako, C. Boscagli, M. Müller, F. Ebert, R. N. d’Alnoncourt, A. Schaefer, S. A. Schunk, I. Jevtovikj, F. Rosowski and S. De, Catal. Sci. Technol., 2023, 13, 2656–2661.

---