Operando Modeling of Zeolite-Catalyzed Reactions Using First-Principles Molecular Dynamics Simulations

Within this Perspective, we critically reflect on the role of first-principles molecular dynamics (MD) simulations in unraveling the catalytic function within zeolites under operating conditions. First-principles MD simulations refer to methods where the dynamics of the nuclei is followed in time by integrating the Newtonian equations of motion on a potential energy surface that is determined by solving the quantum-mechanical many-body problem for the electrons. Catalytic solids used in industrial applications show an intriguing high degree of complexity, with phenomena taking place at a broad range of length and time scales. Additionally, the state and function of a catalyst critically depend on the operating conditions, such as temperature, moisture, presence of water, etc. Herein we show by means of a series of exemplary cases how first-principles MD simulations are instrumental to unravel the catalyst complexity at the molecular scale. Examples show how the nature of reactive species at higher catalytic temperatures may drastically change compared to species at lower temperatures and how the nature of active sites may dynamically change upon exposure to water. To simulate rare events, first-principles MD simulations need to be used in combination with enhanced sampling techniques to efficiently sample low-probability regions of phase space. Using these techniques, it is shown how competitive pathways at operating conditions can be discovered and how broad transition state regions can be explored. Interestingly, such simulations can also be used to study hindered diffusion under operating conditions. The cases shown clearly illustrate how first-principles MD simulations reveal insights into the catalytic function at operating conditions, which could not be discovered using static or local approaches where only a few points are considered on the potential energy surface (PES). Despite these advantages, some major hurdles still exist to fully integrate first-principles MD methods in a standard computational catalytic workflow or to use the output of MD simulations as input for multiple length/time scale methods that aim to bridge to the reactor scale. First of all, methods are needed that allow us to evaluate the interatomic forces with quantum-mechanical accuracy, albeit at a much lower computational cost compared to currently used density functional theory (DFT) methods. The use of DFT limits the currently attainable length/time scales to hundreds of picoseconds and a few nanometers, which are much smaller than realistic catalyst particle dimensions and time scales encountered in the catalysis process. One solution could be to construct machine learning potentials (MLPs), where a numerical potential is derived from underlying quantum-mechanical data, which could be used in subsequent MD simulations. As such, much longer length and time scales could be reached; however, quite some research is still necessary to construct MLPs for the complex systems encountered in industrially used catalysts. Second, most currently used enhanced sampling techniques in catalysis make use of collective variables (CVs), which are mostly determined based on chemical intuition. To explore complex reactive networks with MD simulations, methods are needed that allow the automatic discovery of CVs or methods that do not rely on a priori definition of CVs. Recently, various data-driven methods have been proposed, which could be explored for complex catalytic systems. Lastly, first-principles MD methods are currently mostly used to investigate local reactive events. We hope that with the rise of data-driven methods and more efficient methods to describe the PES, first-principles MD methods will in the future also be able to describe longer length/time scale processes in catalysis. This might lead to a consistent dynamic description of all steps-diffusion, adsorption, and reaction-as they take place at the catalyst particle level.

Constructing Collective Variables Using Invariant Learned Representations

On the time scales accessible to atomistic numerical modeling, chemical reactions are considered rare events. Therefore, the atomistic simulations are commonly biased along a low-dimensional representation of a chemical reaction in an atomic structure space, i.e., along the collective variables. However, suitable collective variables are often complicated to guess a priori. We propose a novel method of collective variable discovery based on dimensionality reduction of the atomic representation vectors. These linear-scaling and invariant representations can be either fixed (untrained) or learned by supervised training of the end-to-end machine learning potential. The learned representations are expected to reflect not only the structural but also the energetic features of the system that are transferable to all of the reactive transformation covered by the machine learning potential. We demonstrate our approach on four high-barrier reactions ranging from a simple gas-phase hydrogen jump reaction to complex reactions in periodic models of industrially relevant heterogeneous catalysts. High data efficiency, automatized feature extraction, favorable scaling, and retention of inherent invariances are all properties that are expected to enable fast and largely automatic construction of suitable collective variables even in highly complex reactive scenarios such as reactive/catalytic transformations at solid-liquid interfaces.

Mechanistic Characterization of Zeolite-Catalyzed Aromatic Electrophilic Substitution at Realistic Operating Conditions

Zeolite-catalyzed benzene ethylation is an important industrial reaction, as it is the first step in the production of styrene for polymer manufacturing. Furthermore, it is a prototypical example of aromatic electrophilic substitution, a key reaction in the synthesis of many bulk and fine chemicals. Despite extensive research, the reaction mechanism and the nature of elusive intermediates at realistic operating conditions is not properly understood. More in detail, the existence of the elusive arenium ion (better known as Wheland complex) formed upon electrophilic attack on the aromatic ring is still a matter of debate. Temperature effects and the presence of protic guest molecules such as water are expected to impact the reaction mechanism and lifetime of the reaction intermediates. Herein, we used enhanced sampling ab initio molecular dynamics simulations to investigate the complete mechanism of benzene ethylation with ethene and ethanol in the H-ZSM-5 zeolite. We show that both the stepwise and concerted mechanisms are active at reaction conditions and that the Wheland intermediate spontaneously appears as a shallow minimum in the free energy surface after the electrophilic attack on the benzene ring. Addition of water enhances the protonation kinetics by about 1 order of magnitude at coverages of one water molecule per Brønsted acidic site. In the fully solvated regime, an overstabilization of the BAS as hydronium ion occurs and the rate enhancement disappears. The obtained results give critical atomistic insights in the role of water to selectively tune the kinetics of protonation reactions in zeolites.

Ab initio molecular dynamics with enhanced sampling in heterogeneous catalysis

Enhanced sampling ab initio simulations enable to study chemical phenomena in catalytic systems including thermal effects & anharmonicity, & collective dynamics describing enthalpic & entropic contributions, which can significantly impact on reaction free energy landscapes.

Revisiting the Structure of Calcined and Hydrated AlPO-11 with DFT-Based Molecular Dynamics Simulations*

Published crystal structures of the AEL-type aluminophosphate AlPO-11 in its calcined form (space group I m a 2 ) show some peculiar features, such as unusually short Al-O and P-O bonds and near-linear Al-O-P angles. Although experimental evidence for the presence of dynamic disorder was presented, the nature of the associated distortions remained unresolved. In this study, ab initio molecular dynamics (AIMD) calculations in the framework of density functional theory (DFT) were employed to study the dynamic behaviour of this zeotype. At 100 K, static local distortions that break the I m a 2 symmetry are present in the time-averaged structures computed from the AIMD trajectories. At 300 and 500 K, the time-averaged structures approach I m a 2 symmetry. Although shortened Al-O and P-O bonds and near-linear Al-O-P angles were found in the average structures, an analysis of radial and angular distribution functions confirmed their absence in the instantaneous structures. This deviation is due to a precession-like motion of some oxygen atoms around the Al-P connection line, which moves their time-averaged positions closer to this line. In hydrated AlPO-11, some of the water molecules are coordinated to framework Al atoms, leading to an octahedral coordination of 1/5 of the Al sites. DFT optimisations and AIMD simulations on partially hydrated models delivered evidence for a preferential adsorption at the Al1 site. No dynamic disorder was observed for the hydrated form.

Solvation and Mobilization of Copper Active Sites in Zeolites by Ammonia: Consequences for the Catalytic Reduction of Nitrogen Oxides

ConspectusCopper-exchanged chabazite (Cu-CHA) zeolites are catalysts used in diesel emissions control for the abatement of nitrogen oxides (NO_x) via selective catalytic reduction (SCR) reactions with ammonia as the reductant. The discovery of these materials in the early 2010s enabled a step-change improvement in diesel emissions aftertreatment technology. Key advantages of Cu-CHA zeolites over prior materials include their effectiveness at the lower temperatures characteristic of diesel exhaust, their durability under high-temperature hydrothermal conditions, and their resistance to poisoning from residual hydrocarbons present in exhaust. Fundamental catalysis research has since uncovered mechanistic and kinetic features that underpin the ability of Cu-CHA to selectively reduce NO_x under strongly oxidizing conditions and to achieve improved NO_x conversion relative to other zeolite frameworks, particularly at low exhaust temperatures and with ammonia instead of other reductants.One critical mechanistic feature is the NH₃ solvation of exchanged Cu ions at low temperatures (<523 K) to create cationic Cu-amine coordination complexes that are ionically tethered to anionic Al framework sites. This ionic tethering confers regulated mobility that facilitates interconversion between mononuclear and binuclear Cu complexes, which is necessary to propagate SCR through a Cu²⁺/Cu⁺ redox cycle during catalytic turnover. This dynamic catalytic mechanism, wherein single and dual metal sites interconvert to mediate different half-reactions of the redox cycle, combines features canonically associated with homogeneous and heterogeneous reaction mechanisms.In this Account, we describe how a unified experimental and theoretical interrogation of Cu-CHA catalysts in operando provided quantitative evidence of regulated Cu ion mobility and its role in the SCR mechanism. This approach relied on new synthetic methods to prepare model Cu-CHA zeolites with varied active-site structures and spatial densities in order to verify that the kinetic and mechanistic models describe the catalytic behavior of a family of materials of diverse composition, and on new computational approaches to capture the active-site structure and dynamics under conditions representative of catalysis. Ex situ interrogation revealed that the Cu structure depends on the conditions for the zeolite synthesis, which influence the framework Al substitution patterns, and that statistical and electronic structure models can enumerate Cu site populations for a known Al distribution. This recognition unifies seemingly disparate spectroscopic observations and inferences regarding Cu ion structure and responses to different external conditions. SCR rates depend strongly on the Cu spatial density and zeolite composition in kinetic regimes where Cu⁺ oxidation with O₂ becomes rate-limiting, as occurs at lower temperatures and under fuel-rich conditions. Transient experiments, ab initio molecular dynamics simulations, and statistical models relate these sensitivities to the mobility constraints imposed by the CHA framework on NH₃-solvated Cu ions, which regulate the pore volume accessible to these ions and their ability to pair and complete the catalytic cycle. This highlights the key characteristics of the CHA framework that enable superior performance under low-temperature SCR reaction conditions.This work illustrates the power of precise control over a catalytic material, simultaneous kinetic and spectroscopic interrogation over a wide range of reaction conditions, and computational strategies tailored to capture those reaction conditions to reveal in microscopic detail the mechanistic features of a complex and widely practiced catalysis. In doing so, it highlights the key role of ion mobility in catalysis and thus potentially a more general phenomenon of reactant solvation and active site mobilization in reactions catalyzed by exchanged metal ions in zeolites.

Proton Acidity and Proton Mobility in ECR-40, a Silicoaluminophosphate that Violates Löwenstein's Rule

The silicoaluminophosphate zeotype ECR-40 contains linkages of AlO4 tetrahedra via a common oxygen atom, thereby violating the famous "Löwenstein's rule". In this work, a combination of static density functional theory (DFT) calculations and DFT-based ab-initio molecular dynamics (AIMD) simulations were employed to study the acidity and mobility of protons associated with such unusual linkages. It was found that the Al-O-Al linkages are preferentially protonated, as deprotonation causes a local accumulation of negative charge. The protons at these linkages possess a somewhat lower Brønsted acidity than those at Si-O-Al links. AIMD simulations for fully hydrated ECR-40 predicted a partial deprotonation of the Al-O-Al linkages, whereas Si-O-Al linkages were fully deprotonated. Frequently, a coordination of water molecules to framework Al atoms was observed in the vicinity of the Al-O-Al links. Hence, these linkages appear prone to break upon dehydration, potentially explaining why Löwenstein's rule is mostly obeyed in materials formed in aqueous media.

Fast room temperature lability of aluminosilicate zeolites

Aluminosilicate zeolites are traditionally used in high-temperature applications at low water vapour pressures where the zeolite framework is generally considered to be stable and static. Increasingly, zeolites are being considered for applications under milder aqueous conditions. However, it has not yet been established how neutral liquid water at mild conditions affects the stability of the zeolite framework. Here, we show that covalent bonds in the zeolite chabazite (CHA) are labile when in contact with neutral liquid water, which leads to partial but fully reversible hydrolysis without framework degradation. We present ab initio calculations that predict novel, energetically viable reaction mechanisms by which Al-O and Si-O bonds rapidly and reversibly break at 300 K. By means of solid-state NMR, we confirm this prediction, demonstrating that isotopic substitution of ¹⁷O in the zeolitic framework occurs at room temperature in less than one hour of contact with enriched water.

While aluminosilicate zeolites are of interest for many applications, the affect of water on zeolite stability in mild aqueous conditions has yet to be established. Here, using ab initio calculations and NMR spectroscopy, the authors show that covalent bonds in the zeolite chabazite are labile when in contact with neutral liquid water.

The effect of water on the validity of Löwenstein's rule

Löwenstein's rule is explained in terms of the level of solvating water inside zeolite pores, along with the formation of Brønsted acidic water clusters derived from framework sites.

The common understanding of zeolite acidity is based on Löwenstein's rule, which states that Al–O–Al aluminium pairs are forbidden in zeolites. This rule is generally accepted to be inviolate in zeolites. However, recent computational research using a 0 K DFT model has suggested that the rule is violated for the acid form of several zeolites under anhydrous conditions [Fletcher et al., Chem. Sci., 8, (2017), 7483]. The effect of water loading on the preferred aluminium distribution in zeolites, however, has so far not been taken into account. In this article, we show by way of ab initio molecular dynamics simulations that Löwenstein's rule is obeyed under high water solvation for acid chabazite (H-CHA) but disobeyed under anhydrous conditions. We find that varying the water loading in the pores leads to dramatic effects on the structure of the active sites and the dynamics of solvation. The solvation of Brønsted protons in the surrounding water was found to be the energetic driving force for the preferred Löwenstein Al distribution and this driving force is absent in non-Löwenstein (Al–O(H)–Al) moieties. The preference for solvated protons further implies that the catalytically active species in zeolites is a protonated water cluster, rather than a framework Brønsted site. Hence, an accurate treatment of the solvation conditions is crucial to capture the behaviour of zeolites and to properly connect simulations to experiments. This work should lead to a change in modelling paradigm for zeolites, from single molecules towards high solvation models where appropriate.

Ethene Dimerization on Zeolite-Hosted Ni Ions: Reversible Mobilization of the Active Site

The active site in ethene oligomerization catalyzed by Ni-zeolites is proposed to be a mobile Ni(II) complex, based on density functional theory-based molecular dynamics (DFT-MD) simulations corroborated by continuous-flow experiments on Ni-SSZ-24 zeolite. The results of the simulations at operating conditions show that ethene molecules reversibly mobilize the active site as they exchange with the zeolite as ligands on Ni during reaction. Microkinetic modeling was conducted on the basis of free-energy profiles derived with DFT-MD for oligomerization on these mobile [(ethene)₂-Ni-alkyl]⁺ species. The model reproduces the experimentally observed high selectivity to dimerization and indicates that the mechanism is consistent with the observed second-order rate dependence on ethene pressure.

Advances in Sustainable Catalysis: A Computational Perspective

The enormous challenge of moving our societies to a more sustainable future offers several exciting opportunities for computational chemists. The first principles approach to "catalysis by design" will enable new and much greener chemical routes to produce vital fuels and fine chemicals. This prospective outlines a wide variety of case studies to underscore how the use of theoretical techniques, from QM/MM to unrestricted DFT and periodic boundary conditions, can be applied to biocatalysis and to both homogeneous and heterogenous catalysts of all sizes and morphologies to provide invaluable insights into the reaction mechanisms they catalyze.

Kinetic Modeling of Catalytic Olefin Cracking and Methanol-to-Olefins (MTO) over Zeolites: A Review

The increasing demand for lower olefins requires new production routes besides steam cracking and fluid catalytic cracking (FCC). Furthermore, less energy consumption, more flexibility in feed and a higher influence on the product distribution are necessary. In this context, catalytic olefin cracking and methanol-to-olefins (MTO) gain in importance. Here, the undesired higher olefins can be catalytically converted and, for methanol, the possibility of a green synthesis route exists. Kinetic modeling of these processes is a helpful tool in understanding the reactivity and finding optimum operating points; however, it is also challenging because reaction networks for hydrocarbon interconversion are rather complex. This review analyzes different deterministic kinetic models published in the literature since 2000. After a presentation of the underlying chemistry and thermodynamics, the models are compared in terms of catalysts, reaction setups and operating conditions. Furthermore, the modeling methodology is shown; both lumped and microkinetic approaches can be found. Despite ZSM-5 being the most widely used catalyst for these processes, other catalysts such as SAPO-34, SAPO-18 and ZSM-23 are also discussed here. Finally, some general as well as reaction-specific recommendations for future work on modeling of complex reaction networks are given.

Towards operando computational modeling in heterogeneous catalysis

An increased synergy between experimental and theoretical investigations in heterogeneous catalysis has become apparent during the last decade. Experimental work has extended from ultra-high vacuum and low temperature towards operando conditions. These developments have motivated the computational community to move from standard descriptive computational models, based on inspection of the potential energy surface at 0 K and low reactant concentrations (0 K/UHV model), to more realistic conditions. The transition from 0 K/UHV to operando models has been backed by significant developments in computer hardware and software over the past few decades. New methodological developments, designed to overcome part of the gap between 0 K/UHV and operando conditions, include (i) global optimization techniques, (ii) ab initio constrained thermodynamics, (iii) biased molecular dynamics, (iv) microkinetic models of reaction networks and (v) machine learning approaches. The importance of the transition is highlighted by discussing how the molecular level picture of catalytic sites and the associated reaction mechanisms changes when the chemical environment, pressure and temperature effects are correctly accounted for in molecular simulations. It is the purpose of this review to discuss each method on an equal footing, and to draw connections between methods, particularly where they may be applied in combination.

Molecular Modelling for Reactor Design

Chemical reactor modelling based on insights and data on a molecular level has become reality over the last few years. Multiscale models describing elementary reaction steps and full microkinetic schemes, pore structures, multicomponent adsorption and diffusion inside pores, and entire reactors have been presented. Quantum mechanical (QM) approaches, molecular simulations (Monte Carlo and molecular dynamics), and continuum equations have been employed for this purpose. Some recent developments in these approaches are presented, in particular time-dependent QM methods, calculation of van der Waals forces, new approaches for force field generation, automatic setup of reaction schemes, and pore modelling. Multiscale simulations are discussed. Applications of these approaches to heterogeneous catalysis are demonstrated for examples that have found growing interest over the last few years, such as metal-support interactions, influence of pore geometry on reactions, noncovalent bonding, reaction dynamics, dynamic changes in catalyst nanoparticle structure, electrocatalysis, solvent effects in catalysis, and multiscale modelling.

On the intrinsic dynamic nature of the rigid UiO-66 metal-organic framework

UiO-66 is a showcase example of an extremely stable metal-organic framework, which maintains its structural integrity during activation processes such as linker exchange and dehydration. The framework can even accommodate a substantial number of defects without compromising its stability. These observations point to an intrinsic dynamic flexibility of the framework, related to changes in the coordination number of the zirconium atoms. Herein we follow the dynamics of the framework in situ, by means of enhanced sampling molecular dynamics simulations such as umbrella sampling, during an activation process, where the coordination number of the bridging hydroxyl groups capped in the inorganic Zr6(μ3-O)4(μ3-OH)4 brick is reduced from three to one. Such a reduction in the coordination number occurs during the dehydration process and in other processes where defects are formed. We observe a remarkable fast response of the system upon structural changes of the hydroxyl group. Internal deformation modes are detected, which point to linker decoordination and recoordination. Detached linkers may be stabilized by hydrogen bonds with hydroxyl groups of the inorganic brick, which gives evidence for an intrinsic dynamic acidity even in the absence of protic guest molecules. Our observations yield a major step forward in the understanding on the molecular level of activation processes realized experimentally but that is hard to track on a purely experimental basis.

Understanding Brønsted-Acid Catalyzed Monomolecular Reactions of Alkanes in Zeolite Pores by Combining Insights from Experiment and Theory

Acidic zeolites are effective catalysts for the cracking of large hydrocarbon molecules into lower molecular weight products required for transportation fuels. However, the ways in which the zeolite structure affects the catalytic activity at Brønsted protons are not fully understood. One way to characterize the influence of the zeolite structure on the catalysis is to study alkane cracking and dehydrogenation at very low conversion, conditions for which the kinetics are well defined. To understand the effects of zeolite structure on the measured rate coefficient (k_app ), it is necessary to identify the equilibrium constant for adsorption into the reactant state (K_ads-H+ ) and the intrinsic rate coefficient of the reaction (k_int ) at reaction temperatures, since k_app is proportional to the product of K_ads-H+ and k_int . We show that K_ads-H+ cannot be calculated from experimental adsorption data collected near ambient temperature, but can, however, be estimated accurately from configurational-bias Monte Carlo (CBMC) simulations. Using monomolecular cracking and dehydrogenation of C₃ -C₆ alkanes as an example, we review recent efforts aimed at elucidating the influence of the acid site location and the zeolite framework structure on the observed values of k_app and its components, K_ads-H+ and k_int .

Metal-organic and covalent organic frameworks as single-site catalysts

Heterogeneous single-site catalysts consist of isolated, well-defined, active sites that are spatially separated in a given solid and, ideally, structurally identical. In this review, the potential of metal-organic frameworks (MOFs) and covalent organic frameworks (COFs) as platforms for the development of heterogeneous single-site catalysts is reviewed thoroughly. In the first part of this article, synthetic strategies and progress in the implementation of such sites in these two classes of materials are discussed. Because these solids are excellent playgrounds to allow a better understanding of catalytic functions, we highlight the most important recent advances in the modelling and spectroscopic characterization of single-site catalysts based on these materials. Finally, we discuss the potential of MOFs as materials in which several single-site catalytic functions can be combined within one framework along with their potential as powerful enzyme-mimicking materials. The review is wrapped up with our personal vision on future research directions.

Towards molecular control of elementary reactions in zeolite catalysis by advanced molecular simulations mimicking operating conditions

Ingeniously selecting zeolite topology and acidity, reaction temperature and guest molecule loading enables tuning the reaction mechanism of zeolite-catalyzed methylation reactions.

Interaction of water with (silico)aluminophosphate zeotypes: a comparative investigation using dispersion-corrected DFT

The adsorption of water in six structurally different aluminophosphates and their silicoaluminophosphate analogues is investigated using dispersion-corrected density-functional theory calculations. In addition to predicting the interaction energies, the structural changes of the materials upon water adsorption are assessed.

Water adsorption in SAPO-34: elucidating the role of local heterogeneities and defects using dispersion-corrected DFT calculations

The chabazite-type silicoaluminophosphate SAPO-34 is a promising adsorbent for applications in thermal energy storage using water adsorption-desorption cycles. In order to develop a microscopic understanding of the impact of local heterogeneities and defects on the water adsorption properties, the interaction of different models of SAPO-34 with water was studied using dispersion-corrected density-functional theory (DFT-D) calculations. In addition to SAPO-34 with isolated silicon atoms, the calculations considered models incorporating two types of heterogeneities (silicon islands, aluminosilicate domains), and two defect-containing (partially and fully desilicated) systems. DFT-D optimisations were performed for systems with small amounts of adsorbed water, in which all H2O molecules can interact with framework protons, and systems with large amounts of adsorbed water (30 H2O molecules per unit cell). At low loadings, the host-guest interaction energy calculated for SAPO-34 with isolated Si atoms amounts to approximately -90 kJ mol(-1). While the presence of local heterogeneities leads to the creation of some adsorption sites that are energetically slightly more favourable, the interaction strength is drastically reduced in systems with defects. At high water loadings, energies in the range of -70 kJ mol(-1) are obtained for all models. The DFT-D interaction energies are in good agreement with experimentally measured heats of water adsorption. A detailed analysis of the equilibrium structures was used to gain insights into the binding modes at low coverages, and to assess the extent of framework deprotonation and changes in the coordination environment of aluminium atoms at high water loadings.